~maddevelopers/mg5amcnlo/simple_unlops

import madgraph

if madgraph.ordering:

    set = misc.OrderedSet

        filename = 'maxconfigs.inc'

        self.write_maxconfigs_file(writers.FortranWriter(filename),

                max(nconfigs,matrix_element.born_me.get_number_of_amplitudes()))

        filename = 'a0Gmuconv.inc'

        startfroma0 = self.write_a0gmuconv_file(

                            writers.FortranWriter(filename),

                            matrix_element)

 

        filename = 'rescale_alpha_tagged.f'

        self.write_rescale_a0gmu_file(

                            writers.FortranWriter(filename),

                            startfroma0, matrix_element)

 

        

        filename = 'maxparticles.inc'

        self.write_maxparticles_file(writers.FortranWriter(filename),

                                     nexternal)

        

    def write_a0gmuconv_file(self, writer, matrix_element):

        """writes an include file with the informations about the 

        alpha0 < > gmu conversion, to be used when the process has

        tagged photons

        """

 

        bool_dict = {True: '.true.', False: '.false.'}

        bornproc = matrix_element.born_me['processes'][0]

        startfromalpha0 = False

        if any([l['is_tagged'] and l['id'] == 22 for l in bornproc['legs']]):

            if 'loop_qcd_qed_sm_a0' in bornproc['model'].get('modelpath'):

                startfromalpha0 = True

 

        text = 'logical  startfroma0\nparameter (startfroma0=%s)\n' % bool_dict[startfromalpha0]

        writer.writelines(text)

        return startfromalpha0

 

    def write_rescale_a0gmu_file(self, writer, startfroma0, matrix_element):

        """writes the function that computes the rescaling factor needed in

        the case of external photons

        """

 

        # get the model parameters

        params = sum([v for v in self.model.get('parameters').values()], [])

        parnames = [p.name.lower() for p in params]

 

        bornproc = matrix_element.born_me['processes'][0]

        # this is to ensure compatibility with standard processes

        if not any([l['is_tagged'] and l['id'] == 22 for l in bornproc['legs']]):

            to_check = []

            expr = '1d0'

            conv_pol = '0d0'

            conv_fin = '0d0'

        

        elif startfroma0:

            to_check = ['mdl_aewgmu', 'mdl_aew']

            base = 'mdl_aewgmu/mdl_aew'

            exp = 'qed_pow/2d0-ntag'

            expr = '(%s)**(%s)' % (base, exp)

            conv_fin = '(qed_pow - ntagph * 2d0) * MDL_ECOUP_DGMUA0_UV_EW_FIN_ * born_wgt'

            conv_pol = '(qed_pow - ntagph * 2d0) * MDL_ECOUP_DGMUA0_UV_EW_1EPS_ * born_wgt'

        else:

            to_check = ['mdl_aew', 'mdl_aew0']

            base = 'mdl_aew0/mdl_aew'

            exp = 'ntag'

            expr = '(%s)**(%s)' % (base, exp)

            conv_fin = '- ntagph * 2d0 * MDL_ECOUP_DGMUA0_UV_EW_FIN_ * born_wgt'

            conv_pol = '- ntagph * 2d0 * MDL_ECOUP_DGMUA0_UV_EW_1EPS_ * born_wgt'

 

        replace_dict = {'rescale_fact': expr,

                        'virtual_a0Gmu_conv_finite': conv_fin,

                        'virtual_a0Gmu_conv_pole': conv_pol}

 

        if not all(p in parnames for p in to_check):

            raise fks_common.FKSProcessError(

                    'Some parameters needed when there are tagged '+\

                    'photons cannot be found in the model.\n' +\

                    'Please load the correct model and restriction ' +\

                    '(e.g loop_qcd_qed_sm_Gmu-a0 or loop_qcd_qed_sm_a0-Gmu)')

 

        file = open(os.path.join(_file_path, \

                          'iolibs/template_files/rescale_alpha_tagged.inc')).read()

        file = file % replace_dict

        

        # Write the file

        writer.writelines(file)

 

 

    def write_maxparticles_file(self, writer, maxparticles):

        """Write the maxparticles.inc file for MadEvent"""

        lines = "integer max_particles, max_branch\n"

        lines += "parameter (max_particles=%d) \n" % maxparticles

        lines += "parameter (max_branch=max_particles-1)"

        writer.writelines(lines)

 

 

    def write_maxconfigs_file(self, writer, maxconfigs):

        """Write the maxconfigs.inc file for MadEvent"""

        lines = "integer lmaxconfigs\n"

        lines += "parameter (lmaxconfigs=%d)" % maxconfigs

        writer.writelines(lines)

 

 

 

        # Add information relevant for FxFx matching:

        # Maximum QCD power in all the contributions

        max_qcd_order = 0

        for diag in matrix_element.get('diagrams'):

            orders = diag.calculate_orders()

            if 'QCD' in orders:

                max_qcd_order = max(max_qcd_order,orders['QCD'])  

        max_n_light_final_partons = max(len([1 for id in proc.get_final_ids() 

        if proc.get('model').get_particle(id).get('mass')=='ZERO' and

               proc.get('model').get_particle(id).get('color')>1])

                                    for proc in matrix_element.get('processes'))

        # Maximum number of final state light jets to be matched

        self.proc_characteristic['max_n_matched_jets'] = max(

                               self.proc_characteristic['max_n_matched_jets'],

                                   min(max_qcd_order,max_n_light_final_partons))   

 

        misc.sprint(self.proc_characteristic['max_n_matched_jets'], max_qcd_order,max_n_light_final_partons)

            tag_lines = []

                tag_lines.append( \

                    'DATA (PARTICLE_TAG_D(%d, IPOS), IPOS=1, NEXTERNAL) / %s /' \

                    % (i + 1, ', '.join(bool_dict[tag] for tag in fksborn.real_processes[info['n_me']-1].particle_tags) ))

            tags = [l.get('is_tagged') for l in bornproc.get('legs')] + [False]

            tag_lines = ['DATA (PARTICLE_TAG_D(1, IPOS), IPOS=1, NEXTERNAL) / %s /' \

                            %  ', '.join(bool_dict[tag] for tag in tags)]

        replace_dict['tag_lines'] = '\n'.join(tag_lines)

        stop 1

        return

        endif

        end

        """

        

        lines+= """

        function get_spin(ipdg)

        implicit none

        integer get_spin, ipdg

 

        if(ipdg.eq.%d)then

        get_spin=%d

        return

        """ % (particle_ids[0], model.get_particle(particle_ids[0]).get('spin'))

 

        for part_id in particle_ids[1:]:

            lines += """else if(ipdg.eq.%d)then

            get_spin=%d

            return

            """ % (part_id, model.get_particle(part_id).get('spin'))

        # Dummy particle for multiparticle vertices with pdg given by

        # first code not in the model

        lines += """else if(ipdg.eq.%d)then

c           This is dummy particle used in multiparticle vertices

            get_spin=-2

            return

            """ % model.get_first_non_pdg()

        lines += """else

        write(*,*)'Error: No spin given for pdg ',ipdg

        stop 1