~mg5core1/mg5amcnlo/2.6.4

« back to all changes in this revision

Viewing changes to madgraph/various/process_checks.py

1. Fixed many of the first series of points from Olivier's review

Show diffs side-by-side

added added

removed removed

Lines of Context:
30
30
import shutil
31
31
import glob
32
32
import re
 
33
import subprocess
33
34
 
34
35
import aloha.aloha_writers as aloha_writers
35
36
import aloha.create_aloha as create_aloha
390
391
            import madgraph.loop.loop_exporters as loop_exporters
391
392
            FortranExporter = loop_exporters.LoopProcessExporterFortranSA(\
392
393
                            self.mg_root, export_dir, True,\
393
 
                            os.path.join(self.mg_root, 'loop_material'),\
 
394
                            os.path.join(self.mg_root, 'Template/loop_material'),\
394
395
                            self.cuttools_dir)
395
396
            FortranModel = helas_call_writers.FortranUFOHelasCallWriter(model)
396
397
            FortranExporter.copy_v4template(modelname=model.get('name'))
402
403
            FortranExporter.finalize_v4_directory(helas_objects.HelasMatrixElementList([\
403
404
                                    matrix_element]),"",False,False,'gfortran')
404
405
 
405
 
        # Same comment as above for the import of LoopMG5Runner
406
 
        from tests.parallel_tests.loop_me_comparator import LoopMG5Runner
407
 
        LoopMG5Runner.fix_PSPoint_in_check(export_dir)
 
406
        self.fix_PSPoint_in_check(export_dir)
408
407
        if gauge_check:
409
408
            file_path, orig_file_content, new_file_content = \
410
409
                                            self.setup_ward_check(export_dir)
412
411
            file.write(new_file_content)
413
412
            file.close()
414
413
        
415
 
        finite_m2 = LoopMG5Runner.get_me_value(process.shell_string_v4(), 0,\
 
414
        finite_m2 = self.get_me_value(process.shell_string_v4(), 0,\
416
415
                                               export_dir, p,verbose=False)[0][0]
417
416
        
418
417
        # Restore the original loop_matrix.f code so that it could be reused
433
432
        else:
434
433
            return {'m2': finite_m2, output:[]}
435
434
 
 
435
    def fix_PSPoint_in_check(self,dir_name):
 
436
        """Set check_sa.f to be reading PS.input assuming a working dir dir_name"""
 
437
 
 
438
        file = open(os.path.join(dir_name, 'SubProcesses', 'check_sa.f'), 'r')
 
439
        check_sa = file.read()
 
440
        file.close()
 
441
 
 
442
        file = open(os.path.join(dir_name, 'SubProcesses', 'check_sa.f'), 'w')
 
443
        file.write(re.sub("READPS = .FALSE.", "READPS = .TRUE.", check_sa))
 
444
        file.close()
 
445
 
 
446
    def get_me_value(self, proc, proc_id, working_dir, PSpoint=[], verbose=True):
 
447
        """Compile and run ./check, then parse the output and return the result
 
448
        for process with id = proc_id and PSpoint if specified."""  
 
449
        if verbose:
 
450
            sys.stdout.write('.')
 
451
            sys.stdout.flush()
 
452
         
 
453
        shell_name = None
 
454
        directories = glob.glob(os.path.join(working_dir, 'SubProcesses',
 
455
                                  'P%i_*' % proc_id))
 
456
        if directories and os.path.isdir(directories[0]):
 
457
            shell_name = os.path.basename(directories[0])
 
458
 
 
459
        # If directory doesn't exist, skip and return 0
 
460
        if not shell_name:
 
461
            logging.info("Directory hasn't been created for process %s" %proc)
 
462
            return ((0.0, 0.0, 0.0, 0.0, 0), [])
 
463
 
 
464
        if verbose: logging.info("Working on process %s in dir %s" % (proc, shell_name))
 
465
        
 
466
        dir_name = os.path.join(working_dir, 'SubProcesses', shell_name)
 
467
        # Run make
 
468
        devnull = open(os.devnull, 'w')
 
469
        retcode = subprocess.call('make',
 
470
                        cwd=dir_name,
 
471
                        stdout=devnull, stderr=devnull)
 
472
                        
 
473
        if retcode != 0:
 
474
            logging.info("Error while executing make in %s" % shell_name)
 
475
            return ((0.0, 0.0, 0.0, 0.0, 0), [])
 
476
 
 
477
        # If a PS point is specified, write out the corresponding PS.input
 
478
        if PSpoint:
 
479
            PSfile = open(os.path.join(dir_name, 'PS.input'), 'w')
 
480
            PSfile.write('\n'.join([' '.join(['%.16E'%pi for pi in p]) \
 
481
                                  for p in PSpoint]))
 
482
            PSfile.close()
 
483
        
 
484
        # Run ./check
 
485
        try:
 
486
            output = subprocess.Popen('./check',
 
487
                        cwd=dir_name,
 
488
                        stdout=subprocess.PIPE, stderr=subprocess.STDOUT).stdout
 
489
            output.read()
 
490
            output.close()
 
491
            if os.path.exists(os.path.join(dir_name,'result.dat')):
 
492
                return self.parse_check_output(file(dir_name+'/result.dat'))  
 
493
            else:
 
494
                logging.warning("Error while looking for file %s"%str(os.path\
 
495
                                                  .join(dir_name,'result.dat')))
 
496
                return ((0.0, 0.0, 0.0, 0.0, 0), [])
 
497
        except IOError:
 
498
            logging.warning("Error while executing ./check in %s" % shell_name)
 
499
            return ((0.0, 0.0, 0.0, 0.0, 0), [])
 
500
 
 
501
    def parse_check_output(self,output):
 
502
        """Parse the output string and return a pair where first four values are 
 
503
        the finite, born, single and double pole of the ME and the fourth is the
 
504
        GeV exponent and the second value is a list of 4 momenta for all particles 
 
505
        involved."""
 
506
 
 
507
        res_p = []
 
508
        value = [0.0,0.0,0.0,0.0]
 
509
        gev_pow = 0
 
510
 
 
511
        for line in output:
 
512
            splitline=line.split()
 
513
            if splitline[0]=='PS':
 
514
                res_p.append([float(s) for s in splitline[1:]])
 
515
            elif splitline[0]=='BORN':
 
516
                value[1]=float(splitline[1])
 
517
            elif splitline[0]=='FIN':
 
518
                value[0]=float(splitline[1])
 
519
            elif splitline[0]=='1EPS':
 
520
                value[2]=float(splitline[1])
 
521
            elif splitline[0]=='2EPS':
 
522
                value[3]=float(splitline[1])
 
523
            elif splitline[0]=='EXP':
 
524
                gev_pow=int(splitline[1])
 
525
 
 
526
        return ((value[0],value[1],value[2],value[3],gev_pow), res_p)
436
527
    
437
528
    def setup_ward_check(self, working_dir):
438
529
        """ Modify loop_matrix.f so to have one external massless gauge boson