#!/usr/bin/python import os,sys,re def steal_numbers(my_string): return re.findall(r'(\d+)', my_string) def ec_r(r): if r<8000 or r>9000: print 'ERROR: r outside of range; r=',r sys.exit(-1) return None def megaparse(fsim,fib,fh71): results = {} # sim log file first for l in file(fsim,'r').readlines(): if l.startswith('#'): continue else: ls = l.split() r = int(ls[0]) #use the GEANT run number as identifier ec_r(r) results[r] = {} for e in ls: if '.dat' in e: results[r]['s'] = e.split('/')[-1] if '00012' in e: field = e.split('/')[-1] if field == 'bfld_map_0_0.00012' or field == 'bfld_map_tuned.00012': results[r]['f'] = 'corrected' elif field == 'bfld_map.00012' or field == 'bfld_map_opera_tuned.00012': results[r]['f'] = 'opera' else: print 'ERROR: do not recognise field entry in this file:' print e sys.exit(-1) # then internal beam table for l in file(fib,'r').readlines(): if l.startswith('#'): continue else: ls = l.split() r = int(ls[0][2:6]) ec_r(r) results[r]['A'] = float(ls[2]) results[r]['L'] = float(ls[4]) results[r]['ph'] = float(ls[6]) results[r]['dx'] = float(ls[7]) results[r]['dy'] = float(ls[8]) # then h71 (output of hm) for l in file(fh71,'r').readlines(): if l.startswith('#'): continue else: ls = l.split() r = int(ls[0]) ec_r(r) results[r]['pmu'] = float(ls[2]) return results def lookup(results,s,f,x,y,tx,ty): req = '_'+str(x)+'_'+str(y)+'_'+str(tx)+'_'+str(ty) for k in results.keys(): z = results[k] if z['s'][:3] == s and req in z['s'] and z['f']==f: return z return None if __name__=='__main__': # read in simulation stuff no_trans = megaparse('/twist/data9/jbueno/simulationOutput/syst_corr_field_no_trans/log.txt', '/twist/data9/jbueno/analysisOutput/syst_corr_field_no_trans/ib_params.txt', '/twist/data9/jbueno/simulationOutput/syst_corr_field_no_trans/h71.txt') inc_trans = megaparse('/twist/data33/jbueno/simulationOutput/syst_corr_field_inc_trans/log.txt', '/twist/data33/jbueno/analysisOutput/syst_corr_field_inc_trans/ib_params.txt', '/twist/data33/jbueno/simulationOutput/syst_corr_field_inc_trans/h71.txt') #debug: # ks = no_trans.keys() # ks.sort() # for k in ks: # print k,no_trans[k] ## __________ following code is for Dx vs Dy plots __________ # s = 's72' # f = 'corrected' # print '#dummy_label TECpos_dx00 TECpos_dy00 TECpos_dx22 TECpos_dy22 dummy_label2 TECang_dx00 TECang_dy00 TECang_dx22 TECang_dy22 ' # for (x,y,tx,ty) in [ (2,0,3,0), (0,-2,0,-3), (-2,0,-3,0), (0,2,0,3), (2,0,3,0) ]: # r1 = lookup(no_trans,s,f,x,y,0,0) # r2 = lookup(inc_trans,s,f,x,y,0,0) # print ('TECpos_'+str(x)+'_'+str(y)).ljust(15),' ',r1['dx'],' ',r1['dy'],' ',r2['dx'],' ',r2['dy'],' ', # r1 = lookup(no_trans,s,f,0,0,tx,ty) # r2 = lookup(inc_trans,s,f,0,0,tx,ty) # print ('TECang_'+str(tx)+'_'+str(ty)).ljust(15),' ',r1['dx'],' ',r1['dy'],' ',r2['dx'],' ',r2['dy'] ## __________ following code is for pmu vs A plots __________ # s = 's72' # f = 'corrected' # print '#dummy_label TECpos00_A TECpos00_pmu TECpos22_A TECpos22_pmu dummy_label2 TECang00_A TECang00_pmu TECang22_A TECang22_pmu ' # for (x,y,tx,ty) in [ (2,0,3,0), (0,-2,0,-3), (-2,0,-3,0), (0,2,0,3)]: #, (2,0,3,0) ]: # r1 = lookup(no_trans,s,f,x,y,0,0) # r2 = lookup(inc_trans,s,f,x,y,0,0) # print ('TECpos_'+str(x)+'_'+str(y)).ljust(15),' ',r1['A'],' ',r1['pmu'],' ',r2['A'],' ',r2['pmu'],' ', # r1 = lookup(no_trans,s,f,0,0,tx,ty) # r2 = lookup(inc_trans,s,f,0,0,tx,ty) # print ('TECang_'+str(tx)+'_'+str(ty)).ljust(15),' ',r1['A'],' ',r1['pmu'],' ',r2['A'],' ',r2['pmu'] ## __________ following code is for pmu vs phi plots __________ # s = 's76' # f = 'opera' # print '#dummy_label TECpos00_ph TECpos00_pmu TECpos22_ph TECpos22_pmu dummy_label2 TECang00_ph TECang00_pmu TECang22_ph TECang22_pmu ' # for (x,y,tx,ty) in [ (2,0,3,0), (0,-2,0,-3), (-2,0,-3,0), (0,2,0,3)]: #, (2,0,3,0) ]: # r1 = lookup(no_trans,s,f,x,y,0,0) # r2 = lookup(inc_trans,s,f,x,y,0,0) # print ('TECpos_'+str(x)+'_'+str(y)).ljust(15),' ',r1['ph'],' ',r1['pmu'],' ',r2['ph'],' ',r2['pmu'],' ', # r1 = lookup(no_trans,s,f,0,0,tx,ty) # r2 = lookup(inc_trans,s,f,0,0,tx,ty) # print ('TECang_'+str(tx)+'_'+str(ty)).ljust(15),' ',r1['ph'],' ',r1['pmu'],' ',r2['ph'],' ',r2['pmu'] ## __________ following code is for pmu vs phi plots __________ # s = 's74' # print '#dummy_label TECpos00_L TECpos00_pmu TECpos22_L TECpos22_pmu dummy_label2 TECang00_L TECang00_pmu TECang22_L TECang22_pmu ' # for (x,y,tx,ty) in [ (2,0,3,0), (0,-2,0,-3), (-2,0,-3,0), (0,2,0,3)]: #, (2,0,3,0) ]: # r1 = lookup(no_trans,s,x,y,0,0) # r2 = lookup(inc_trans,s,x,y,0,0) # print ('TECpos_'+str(x)+'_'+str(y)).ljust(15),' ',r1['L'],' ',r1['pmu'],' ',r2['L'],' ',r2['pmu'],' ', # r1 = lookup(no_trans,s,0,0,tx,ty) # r2 = lookup(inc_trans,s,0,0,tx,ty) # print ('TECang_'+str(tx)+'_'+str(ty)).ljust(15),' ',r1['L'],' ',r1['pmu'],' ',r2['L'],' ',r2['pmu'] ## __________ following code is for pmu vs A plots, independent of what caused it __________ # s = 's86' # f = 'corrected' # print '#A pmu0' # for (x,y,tx,ty) in [ (2,0,3,0), (0,-2,0,-3), (-2,0,-3,0), (0,2,0,3)]: #, (2,0,3,0) ]: # r1 = lookup(no_trans,s,f,x,y,0,0) # r2 = lookup(inc_trans,s,f,x,y,0,0) # print r1['A'],' ',r1['pmu'],' ' # print r2['A'],' ',r2['pmu'],' ' # r1 = lookup(no_trans,s,f,0,0,tx,ty) # r2 = lookup(inc_trans,s,f,0,0,tx,ty) # print r1['A'],' ',r1['pmu'],' ' # print r2['A'],' ',r2['pmu'] # ## __________ following code is for pmu vs A plots, independent of what caused it __________ # s = 's72' #special because there was no tec calib # print '#A pmu0' # for (x,y,tx,ty) in [ (2,0,3,0), (0,-2,0,-3), (-2,0,-3,0), (0,2,0,3)]: #, (2,0,3,0) ]: ## for (x,y,tx,ty) in [ (2,0,6,0), (0,-2,0,-6), (-2,0,-6,0), (0,2,0,6)]: #, (2,0,3,0) ]: # r1 = lookup(no_trans,s,x,y,0,0) # r2 = lookup(inc_trans,s,x,y,0,0) # print r1['A'],' ',r1['pmu'],' ' # print r2['A'],' ',r2['pmu'],' ' # r1 = lookup(no_trans,s,0,0,tx,ty) # r2 = lookup(inc_trans,s,0,0,tx,ty) # print r1['A'],' ',r1['pmu'],' ' # print r2['A'],' ',r2['pmu'] ## __________ following is get dependence of pmu on pos or ang __________ s = 's72' f = 'opera' # f = 'corrected' #following code gets before/after trans # tr = inc_trans # with_trans = lookup(tr,s,f,0,0,0,0)['pmu'] # print with_trans # tr = no_trans # without_trans = lookup(tr,s,f,0,0,0,0)['pmu'] # print without_trans # print 1E4*(with_trans-without_trans) # sys.exit(-1) # poss = [-4,-2,0,2,4] poss = [-2,0,2,] # angs = [-6,-3,0,3,6] angs = [-3,0,3] all_pmu0s = [] tr = inc_trans # tr = no_trans print '#x xpmu y ypmu tx txpmu ty typmu' for i in range(3): x, y, tx, ty = poss[i], poss[i], angs[i], angs[i] r = lookup(tr,s,f,x,0,0,0) print x/10.,r['pmu'],' ', all_pmu0s.append(r['pmu']) r = lookup(tr,s,f,0,y,0,0) print y/10.,r['pmu'],' ', all_pmu0s.append(r['pmu']) r = lookup(tr,s,f,0,0,tx,0) print tx,r['pmu'],' ', all_pmu0s.append(r['pmu']) r = lookup(tr,s,f,0,0,0,ty) print ty,r['pmu'] all_pmu0s.append(r['pmu']) central_value = lookup(tr,s,f,0,0,0,0)['pmu'] min_val = min(all_pmu0s) max_val = max(all_pmu0s) print '\n','-',1E4*(central_value-min_val),' ','+',1E4*(max_val-central_value) # print all_pmu0s # print min(all_pmu0s) # print max(all_pmu0s) # print 1E4*0.5*(min(all_pmu0s)-max(all_pmu0s)) # read in data # data = {'s74':{'dx':0.22, 'dy':0.015}, # 's75':{'dx':0.265,'dy':0.085}, # 's83':{'dx':0.255,'dy':0.04}, # 's84':{'dx':0.265,'dy':-0.05}, # 's87':{'dx':0.245,'dy':0.125}} # results_2006 = megaparse('/twist/data9/jbueno/simulationOutput/bf_trans_brute_2006/log.txt', # '/twist/data9/jbueno/analysisOutput/bf_trans_brute_2006/ib_params.txt') # results_2007 = megaparse('/twist/data33/jbueno/simulationOutput/bf_trans_brute_2007/log.txt', # '/twist/data33/jbueno/analysisOutput/bf_trans_brute_2007/ib_params.txt') # dx,dy = fetch(results_2007,bx,by,'msor6_production_s87_end_set.dat')