ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/RateMonShiftTool_dev/DatabaseRatePredictor.py
Revision: 1.16
Committed: Mon Mar 19 16:13:02 2012 UTC (13 years, 1 month ago) by grchrist
Content type: text/x-python
Branch: MAIN
Changes since 1.15: +26 -24 lines
Log Message:
still messing with RatePredictor, but added StreamA ceiling in cfg file

File Contents

# User Rev Content
1 grchrist 1.3 #!/usr/bin/env python
2    
3 abrinke1 1.1 from DatabaseParser import *
4 abrinke1 1.5 from GetListOfRuns import *
5 abrinke1 1.1 import sys
6     import os
7     from numpy import *
8     import pickle
9    
10     from ROOT import gROOT, TCanvas, TF1, TGraph, TGraphErrors, TPaveStats, gPad, gStyle
11     from ROOT import TFile, TPaveText
12     from ROOT import gBenchmark
13     import array
14     import math
15 grchrist 1.8 from ReadConfig import RateMonConfig
16 abrinke1 1.1
17 abrinke1 1.5 from selectionParser import selectionParser
18    
19 abrinke1 1.1 def main():
20    
21 abrinke1 1.5 ## Can use any combination of LowestRunNumber, HighestRunNumber, and NumberOfRuns -
22     ## just modify "ExistingRuns.sort" and for "run in ExistingRuns" accordingly
23    
24     LowestRunNumber = 176000
25     HighestRunNumber = 180252
26     NumberOfRuns = 1
27    
28     run_list = []
29     ExistingRuns = GetLatestRunNumber(LowestRunNumber)
30     #ExistingRuns.sort(reverse = True) ##Allows you to count down from last run
31    
32     for run in ExistingRuns:
33     #if NumberOfRuns > 0:
34     if run <= HighestRunNumber:
35     run_list.append(run)
36     NumberOfRuns-=1
37    
38 grchrist 1.8 ##
39     ## to get list of triggers
40     ##
41     Config = RateMonConfig(os.path.abspath(os.path.dirname(sys.argv[0])))
42     Config.CFGfile="defaults.cfg"
43     Config.ReadCFG()
44    
45 abrinke1 1.7 ## ###### TO CREATE FITS #########
46 grchrist 1.16 run_list = [179497,179547,179558,179563,179889,179959,179977,180072,180076,180093,180241,180250,180252]
47     #run_list = [180250]
48 grchrist 1.8 trig_name = "HLT"
49 grchrist 1.16 trig_list=["HLT_Ele65_CaloIdVT_TrkIdT_v6", "HLT_HT650_v4", "HLT_IsoMu15_eta2p1_LooseIsoPFTau20_v6", "HLT_IsoMu30_eta2p1_v7", "HLT_Jet370_v10", "HLT_MET200_v7", "HLT_Mu8_Ele17_CaloIdT_CaloIsoVL_v8", "HLT_PFMHT150_v17", "HLT_Photon135_v2", "HLT_Photon26_R9IdT_Photon18_CaloIdXL_IsoXL_Mass60_v4"]
50 grchrist 1.14 ##trig_list=Config.MonitorList
51 grchrist 1.16 ##trig_list=["HLT_IsoMu30_eta2p1_v7"]
52     num_ls = 10
53 grchrist 1.15 physics_active_psi = True ##Requires that physics and active be on, and that the prescale column is not 0
54     JSON = [] ##To not use a JSON file, just leave the array empty
55 grchrist 1.16 JSON = GetJSON("Cert_160404-180252_7TeV_PromptReco_Collisions11_JSON.txt") ##Returns array JSON[runs][ls_list]
56 grchrist 1.15
57 abrinke1 1.1 debug_print = False
58 grchrist 1.14 no_versions=False
59 grchrist 1.15 min_rate = 0.1
60 abrinke1 1.1 print_table = False
61 grchrist 1.15 data_clean = True ##Gets rid of anomalous rate points, reqires physics_active_psi (PAP) and deadtime < 20%
62 abrinke1 1.5 ##plot_properties = [varX, varY, do_fit, save_root, save_png, fit_file]
63 grchrist 1.16 plot_properties = [["delivered", "rate", True, True, False, ""]]
64 grchrist 1.15
65 abrinke1 1.1 masked_triggers = ["AlCa_", "DST_", "HLT_L1", "HLT_L2", "HLT_Zero"]
66 grchrist 1.15 save_fits = True
67 grchrist 1.16 max_dt=0.08 ## no deadtime cutuse 2.0
68 grchrist 1.11 force_new=True
69     print_info=True
70 grchrist 1.15 SubSystemOff={'All':True,'Mu':False,'HCal':True,'ECal':False,'Tracker':False,'EndCap':False,'Beam':True}
71    
72    
73     ###### TO SEE RATE VS PREDICTION ########
74     ## run_list = [180250]
75     ## trig_name = "HLT"
76     ## trig_list = ["HLT_IsoMu30_eta2p1_v7"]
77     ## ##trig_list=["HLT_Ele65_CaloIdVT_TrkIdT_v6", "HLT_HT650_v4", "HLT_IsoMu15_eta2p1_LooseIsoPFTau20_v6", "HLT_IsoMu30_eta2p1_v7", "HLT_Jet370_v10", "HLT_MET200_v7", "HLT_Mu8_Ele17_CaloIdT_CaloIsoVL_v8", "HLT_PFMHT150_v17", "HLT_Photon135_v2", "HLT_Photon26_R9IdT_Photon18_CaloIdXL_IsoXL_Mass60_v4"]
78     ## ##trig_list=Config.MonitorList
79     ## num_ls = 1
80     ## physics_active_psi = True
81     ## JSON = []
82     ## debug_print = False
83     ## no_versions=False
84     ## min_rate = 1.0
85     ## print_table = False
86 grchrist 1.16 ## data_clean = False
87 grchrist 1.15 ## ##plot_properties = [varX, varY, do_fit, save_root, save_png, fit_file]
88     ## ##plot_properties = [["ls", "rawrate", False, True, False, "Fits/2011/Fit_HLT_10LS_Run176023to180252.pkl"]]
89     ## plot_properties = [["ls", "rawrate", False, True, False, "Fits/2011/Fit_HLT_10LS_Run179497to180252.pkl"]]
90 grchrist 1.12
91 grchrist 1.15 ## masked_triggers = ["AlCa_", "DST_", "HLT_L1", "HLT_L2", "HLT_Zero"]
92     ## save_fits = False
93     ## max_dt=2.0 ## no deadtime cut=2.0
94     ## force_new=True
95     ## print_info=True
96     ## SubSystemOff={'All':True,'Mu':False,'HCal':False,'ECal':False,'Tracker':False,'EndCap':False,'Beam':True}
97 grchrist 1.12
98 grchrist 1.11
99 grchrist 1.12 ## print SubSystemOff.keys()
100     ## print SubSystemOff.values()
101 abrinke1 1.1
102 abrinke1 1.2
103 abrinke1 1.5 ######## END PARAMETERS - CALL FUNCTIONS ##########
104 grchrist 1.12 [Rates,LumiPageInfo]= GetDBRates(run_list, trig_name, trig_list, num_ls, max_dt, physics_active_psi, JSON, debug_print, force_new)
105 grchrist 1.10 ##if not checkLS(Rates,LumiPageInfo,trig_list):
106     ## print "Missing LS!"
107 grchrist 1.13
108    
109     ##for iterator in range(len(Rates["HLT_IsoMu30_eta2p1"]["rawrate"])):
110     ## print iterator, "ls=",Rates["HLT_IsoMu30_eta2p1"]["ls"][iterator],"rate=",round(Rates["HLT_IsoMu30_eta2p1"]["rawrate"][iterator],2)
111 grchrist 1.10
112 grchrist 1.11 MakePlots(Rates, LumiPageInfo, run_list, trig_name, trig_list, num_ls, min_rate, max_dt, print_table, data_clean, plot_properties, masked_triggers, save_fits, debug_print,SubSystemOff, print_info)
113 abrinke1 1.1
114    
115 grchrist 1.12 def GetDBRates(run_list,trig_name,trig_list, num_ls, max_dt, physics_active_psi,JSON,debug_print, force_new):
116 abrinke1 1.1
117     Rates = {}
118 grchrist 1.10 LumiPageInfo={}
119 abrinke1 1.5 ## Save in RefRuns with name dependent on trig_name, num_ls, JSON, and physics_active_psi
120     if JSON:
121 grchrist 1.14 print "Using JSON file"
122 abrinke1 1.5 if physics_active_psi:
123     RefRunNameTemplate = "RefRuns/2011/Rates_%s_%sLS_JPAP.pkl"
124     else:
125     RefRunNameTemplate = "RefRuns/2011/Rates_%s_%sLS_JSON.pkl"
126     else:
127 grchrist 1.14 print "Using Physics and Active ==1"
128 abrinke1 1.5 if physics_active_psi:
129     RefRunNameTemplate = "RefRuns/2011/Rates_%s_%sLS_PAP.pkl"
130     else:
131     RefRunNameTemplate = "RefRuns/2011/Rates_%s_%sLS.pkl"
132 grchrist 1.10
133    
134 abrinke1 1.1 RefRunFile = RefRunNameTemplate % (trig_name,num_ls)
135 abrinke1 1.5 RefRunFileHLT = RefRunNameTemplate % ("HLT",num_ls)
136 grchrist 1.10
137     print "RefRun=",RefRunFile
138     print "RefRunFileHLT",RefRunFileHLT
139 grchrist 1.11 if not force_new:
140     try: ##Open an existing RefRun file with the same parameters and trigger name
141     pkl_file = open(RefRunFile, 'rb')
142     Rates = pickle.load(pkl_file)
143     pkl_file.close()
144     os.remove(RefRunFile)
145     print "using",RefRunFile
146    
147    
148 abrinke1 1.5 except:
149 grchrist 1.11 try: ##Open an existing RefRun file with the same parameters and HLT for trigger name
150     pkl_file = open(RefRunFileHLT)
151     HLTRates = pickle.load(pkl_file)
152     for key in HLTRates:
153     if trig_name in str(key):
154     Rates[key] = HLTRates[key]
155     #print str(RefRunFile)+" does not exist. Creating ..."
156     except:
157     print str(RefRunFile)+" does not exist. Creating ..."
158 grchrist 1.10
159     ## try the lumis file
160     RefLumiNameTemplate = "RefRuns/2011/Lumis_%s_%sLS.pkl"
161 grchrist 1.11 RefLumiFile= RefLumiNameTemplate % ("HLT",num_ls)
162     if not force_new:
163     try:
164     pkl_lumi_file = open(RefLumiFile, 'rb')
165     LumiPageInfo = pickle.load(pkl_lumi_file)
166     pkl_lumi_file.close()
167     os.remove(RefLumiFile)
168     print "using",RefLumiFile
169     except:
170     print str(RefLumiFile)+" doesn't exist. Make it..."
171    
172 abrinke1 1.5 for RefRunNum in run_list:
173 grchrist 1.11
174 abrinke1 1.5 if JSON:
175     if not RefRunNum in JSON:
176     continue
177 abrinke1 1.1 try:
178 grchrist 1.10
179 abrinke1 1.1 ExistsAlready = False
180     for key in Rates:
181     if RefRunNum in Rates[key]["run"]:
182     ExistsAlready = True
183     break
184 grchrist 1.10
185    
186     LumiExistsLAready=False
187     for v in LumiPageInfo.itervalues():
188     #print RefRunNum, v["Run"]
189    
190     if RefRunNum == v["Run"]:
191     LumiExistsAlready=True
192     break
193     if ExistsAlready and LumiExistsAlready:
194 abrinke1 1.1 continue
195 grchrist 1.10
196    
197 abrinke1 1.1 except:
198     print "Getting info for run "+str(RefRunNum)
199    
200     if RefRunNum < 1:
201     continue
202 grchrist 1.10 print "calculating rates and green lumis"
203    
204 abrinke1 1.5 if True: ##Placeholder
205 abrinke1 1.1 if True: #May replace with "try" - for now it's good to know when problems happen
206     RefParser = DatabaseParser()
207     RefParser.RunNumber = RefRunNum
208     RefParser.ParseRunSetup()
209 abrinke1 1.5 RefLumiRangePhysicsActive = RefParser.GetLSRange(1,9999) ##Gets array of all LS with physics and active on
210     RefLumiArray = RefParser.GetLumiInfo() ##Gets array of all existing LS and their lumi info
211 abrinke1 1.2 RefLumiRange = []
212 grchrist 1.9 RefMoreLumiArray = RefParser.GetMoreLumiInfo()#dict with keys as bits from lumisections WBM page and values are dicts with key=LS:value=bit
213    
214 grchrist 1.15
215 abrinke1 1.5 for iterator in RefLumiArray[0]: ##Makes array of LS with proper PAP and JSON properties
216 grchrist 1.11 ##cheap way of getting PSCol None-->0
217 grchrist 1.15 if RefLumiArray[0][iterator]!=1 and RefLumiArray[0][iterator]!=2 and RefLumiArray[0][iterator]!=3 and RefLumiArray[0][iterator]!=4 and RefLumiArray[0][iterator]!=5 and RefLumiArray[0][iterator]!=6 and RefLumiArray[0][iterator]!=7 and RefLumiArray[0][iterator]!=8:
218 grchrist 1.11 RefLumiArray[0][iterator]=0
219 grchrist 1.15
220 grchrist 1.11
221     if not physics_active_psi or (RefLumiArray[5][iterator] == 1 and RefLumiArray[6][iterator] == 1):
222 abrinke1 1.5 if not JSON or RefRunNum in JSON:
223     if not JSON or iterator in JSON[RefRunNum]:
224     RefLumiRange.append(iterator)
225 grchrist 1.15 #print iterator, RefLumiArray[0][iterator], "active=",RefLumiArray[5][iterator],"physics=",RefLumiArray[6][iterator]
226 grchrist 1.9 #print "hbhea=",RefMoreLumiArray['hbhea'][iterator]
227    
228 abrinke1 1.5 try:
229     nls = RefLumiRange[0]
230     LSRange = {}
231     except:
232     print "Run "+str(RefRunNum)+" has no good LS"
233     continue
234     if num_ls > len(RefLumiRange):
235 abrinke1 1.1 print "Run "+str(RefRunNum)+" is too short: from "+str(nls)+" to "+str(RefLumiRange[-1])+", while num_ls = "+str(num_ls)
236     continue
237     while nls < RefLumiRange[-1]-num_ls:
238 abrinke1 1.5 LSRange[nls] = []
239     counter = 0
240     for iterator in RefLumiRange:
241     if iterator >= nls and counter < num_ls:
242     LSRange[nls].append(iterator)
243     counter += 1
244 abrinke1 1.1 nls = LSRange[nls][-1]+1
245 abrinke1 1.5
246 abrinke1 1.1 print "Run "+str(RefRunNum)+" contains LS from "+str(min(LSRange))+" to "+str(max(LSRange))
247 grchrist 1.8 for nls in sorted(LSRange.iterkeys()):
248 grchrist 1.9
249 abrinke1 1.1 TriggerRates = RefParser.GetHLTRates(LSRange[nls])
250 grchrist 1.16 #print nls, RefLumiArray[1][nls], RefLumiArray[2][nls]
251 abrinke1 1.5
252     [inst, live, delivered, dead, pscols] = RefParser.GetAvLumiInfo(LSRange[nls])
253    
254 abrinke1 1.2 physics = 1
255     active = 1
256 abrinke1 1.5 psi = 99
257     for iterator in LSRange[nls]: ##Gets lowest value of physics, active, and psi in the set of lumisections
258 abrinke1 1.2 if RefLumiArray[5][iterator] == 0:
259     physics = 0
260     if RefLumiArray[6][iterator] == 0:
261     active = 0
262 abrinke1 1.5 if RefLumiArray[0][iterator] < psi:
263     psi = RefLumiArray[0][iterator]
264 abrinke1 1.2
265 grchrist 1.10 MoreLumiMulti=LumiRangeGreens(RefMoreLumiArray,LSRange,nls,RefRunNum)
266    
267     #print MoreLumiMulti.keys()
268     # print "\n\n\n"
269    
270    
271 grchrist 1.9
272 abrinke1 1.5 if inst < 0 or live < 0 or delivered < 0:
273     print "Run "+str(RefRunNum)+" LS "+str(nls)+" inst lumi = "+str(inst)+" live lumi = "+str(live)+", delivered = "+str(delivered)+", physics = "+str(physics)+", active = "+str(active)
274 abrinke1 1.2
275 grchrist 1.10
276    
277     LumiPageInfo[nls]=MoreLumiMulti
278     ##print LumiPageInfo[nls]
279     ## try:
280     ## LumiPageInfo[keys].append(values)
281     ## except:
282     ## print "Failed",RefRunNum, nls, keys
283    
284 abrinke1 1.1 for key in TriggerRates:
285 grchrist 1.12 ##if not key in trig_list:
286     ## continue
287    
288     ##if not trig_name in key:
289     ## continue
290 abrinke1 1.1 name = key
291 grchrist 1.12
292 grchrist 1.14 ##if re.match('.*_v[0-9]+',name): ##Removes _v#
293     ## name = name[:name.rfind('_')]
294 grchrist 1.12 if not name in trig_list:
295     continue
296     #print "trigger=",name, trig_list
297    
298 abrinke1 1.1 if not Rates.has_key(name):
299     Rates[name] = {}
300     Rates[name]["run"] = []
301     Rates[name]["ls"] = []
302     Rates[name]["ps"] = []
303     Rates[name]["inst_lumi"] = []
304     Rates[name]["live_lumi"] = []
305     Rates[name]["delivered_lumi"] = []
306     Rates[name]["deadtime"] = []
307     Rates[name]["rawrate"] = []
308     Rates[name]["rate"] = []
309     Rates[name]["rawxsec"] = []
310     Rates[name]["xsec"] = []
311 abrinke1 1.2 Rates[name]["physics"] = []
312     Rates[name]["active"] = []
313 abrinke1 1.5 Rates[name]["psi"] = []
314 grchrist 1.9
315 grchrist 1.10 #for keys, values in MoreLumiMulti.iteritems():
316     # Rates[name][keys] = []
317    
318 grchrist 1.9
319 abrinke1 1.1 [avps, ps, rate, psrate] = TriggerRates[key]
320 grchrist 1.11 #print "TriggerRates=",TriggerRates[key], "key=",key
321 abrinke1 1.1 Rates[name]["run"].append(RefRunNum)
322     Rates[name]["ls"].append(nls)
323     Rates[name]["ps"].append(ps)
324     Rates[name]["inst_lumi"].append(inst)
325     Rates[name]["live_lumi"].append(live)
326     Rates[name]["delivered_lumi"].append(delivered)
327     Rates[name]["deadtime"].append(dead)
328     Rates[name]["rawrate"].append(rate)
329 abrinke1 1.2 if live == 0:
330 abrinke1 1.5 Rates[name]["rate"].append(0)
331 abrinke1 1.2 Rates[name]["rawxsec"].append(0.0)
332     Rates[name]["xsec"].append(0.0)
333     else:
334 abrinke1 1.5 Rates[name]["rate"].append(psrate/(1.0-dead))
335 abrinke1 1.2 Rates[name]["rawxsec"].append(rate/live)
336     Rates[name]["xsec"].append(psrate/live)
337     Rates[name]["physics"].append(physics)
338     Rates[name]["active"].append(active)
339 abrinke1 1.5 Rates[name]["psi"].append(psi)
340 grchrist 1.16 #print iterator, "LS=", nls, "dt=",round(dead,2), "deliv=",delivered, "live=",live
341 grchrist 1.9
342 grchrist 1.10 #for keys, values in MoreLumiMulti.iteritems():
343     # Rates[name][keys].append(values)
344 grchrist 1.9 #print nls, name, keys, values
345     #print " "
346 abrinke1 1.5 #except: ##If we replace "if True:" with "try:"
347 abrinke1 1.1 #print "Failed to parse run "+str(RefRunNum)
348    
349 grchrist 1.10 RateOutput = open(RefRunFile, 'wb') ##Save new Rates[] to RefRuns
350     pickle.dump(Rates, RateOutput, 2)
351     RateOutput.close()
352     LumiOutput = open(RefLumiFile,'wb')
353     pickle.dump(LumiPageInfo,LumiOutput, 2)
354     LumiOutput.close()
355    
356    
357     return [Rates,LumiPageInfo]
358 abrinke1 1.1
359 grchrist 1.11 def MakePlots(Rates, LumiPageInfo, run_list, trig_name, trig_list, num_ls, min_rate, max_dt, print_table, data_clean, plot_properties, masked_triggers, save_fits, debug_print, SubSystemOff, print_info):
360 abrinke1 1.1 min_run = min(run_list)
361     max_run = max(run_list)
362 grchrist 1.11
363 abrinke1 1.5 InputFit = {}
364     OutputFit = {}
365 grchrist 1.14 first_trigger=True
366 abrinke1 1.1
367     RootNameTemplate = "%s_%sLS_Run%sto%s.root"
368     RootFile = RootNameTemplate % (trig_name, num_ls, min_run, max_run)
369    
370 abrinke1 1.5 for varX, varY, do_fit, save_root, save_png, fit_file in plot_properties:
371     if fit_file:
372 abrinke1 1.1 pkl_file = open(fit_file, 'rb')
373 abrinke1 1.5 InputFit = pickle.load(pkl_file)
374 abrinke1 1.1 pkl_file.close()
375     if save_root:
376     try:
377     os.remove(RootFile)
378     except:
379     break
380 grchrist 1.15
381 abrinke1 1.1 for print_trigger in Rates:
382 abrinke1 1.7 ##Limits Rates[] to runs in run_list
383     NewTrigger = {}
384 grchrist 1.8 if not print_trigger in trig_list:
385     continue
386 abrinke1 1.7 for key in Rates[print_trigger]:
387     NewTrigger[key] = []
388     for iterator in range (len(Rates[print_trigger]["run"])):
389     if Rates[print_trigger]["run"][iterator] in run_list:
390     for key in Rates[print_trigger]:
391     NewTrigger[key].append(Rates[print_trigger][key][iterator])
392     Rates[print_trigger] = NewTrigger
393    
394 abrinke1 1.1 meanrawrate = sum(Rates[print_trigger]["rawrate"])/len(Rates[print_trigger]["rawrate"])
395     if not trig_name in print_trigger:
396     continue
397     if meanrawrate < min_rate:
398     continue
399     masked_trig = False
400     for mask in masked_triggers:
401     if str(mask) in print_trigger:
402     masked_trig = True
403     if masked_trig:
404     continue
405    
406 abrinke1 1.5 OutputFit[print_trigger] = {}
407 abrinke1 1.1
408     lowlumi = 0
409 abrinke1 1.7 meanlumi_init = median(Rates[print_trigger]["live_lumi"])
410 abrinke1 1.1 meanlumi = 0
411     highlumi = 0
412     lowxsec = 0
413     meanxsec = 0
414     highxsec = 0
415     nlow = 0
416     nhigh = 0
417 grchrist 1.15
418 abrinke1 1.1 for iterator in range(len(Rates[print_trigger]["rate"])):
419     if Rates[print_trigger]["live_lumi"][iterator] <= meanlumi_init:
420 grchrist 1.11 if ( Rates[print_trigger]["rawrate"][iterator] > 0.04 and Rates[print_trigger]["physics"][iterator] == 1 and Rates[print_trigger]["active"][iterator] == 1 and Rates[print_trigger]["deadtime"][iterator] < max_dt and Rates[print_trigger]["psi"][iterator] > 0 and Rates[print_trigger]["live_lumi"] > 500):
421 abrinke1 1.1 meanxsec+=Rates[print_trigger]["xsec"][iterator]
422     lowxsec+=Rates[print_trigger]["xsec"][iterator]
423     meanlumi+=Rates[print_trigger]["live_lumi"][iterator]
424     lowlumi+=Rates[print_trigger]["live_lumi"][iterator]
425     nlow+=1
426     if Rates[print_trigger]["live_lumi"][iterator] > meanlumi_init:
427 grchrist 1.11 if ( Rates[print_trigger]["rawrate"][iterator] > 0.04 and Rates[print_trigger]["physics"][iterator] == 1 and Rates[print_trigger]["active"][iterator] == 1 and Rates[print_trigger]["deadtime"][iterator] < max_dt and Rates[print_trigger]["psi"][iterator] > 0 and Rates[print_trigger]["live_lumi"] > 500):
428 abrinke1 1.1 meanxsec+=Rates[print_trigger]["xsec"][iterator]
429     highxsec+=Rates[print_trigger]["xsec"][iterator]
430     meanlumi+=Rates[print_trigger]["live_lumi"][iterator]
431     highlumi+=Rates[print_trigger]["live_lumi"][iterator]
432     nhigh+=1
433 abrinke1 1.7 try:
434     meanxsec = meanxsec/(nlow+nhigh)
435     meanlumi = meanlumi/(nlow+nhigh)
436     slopexsec = ( (highxsec/nhigh) - (lowxsec/nlow) ) / ( (highlumi/nhigh) - (lowlumi/nlow) )
437     except:
438     print str(print_trigger)+" has no good datapoints - setting initial xsec slope estimate to 0"
439     meanxsec = median(Rates[print_trigger]["xsec"])
440     meanlumi = median(Rates[print_trigger]["live_lumi"])
441     slopexsec = 0
442 abrinke1 1.1
443 abrinke1 1.5 [run_t,ls_t,ps_t,inst_t,live_t,delivered_t,deadtime_t,rawrate_t,rate_t,rawxsec_t,xsec_t,psi_t,e_run_t,e_ls_t,e_ps_t,e_inst_t,e_live_t,e_delivered_t,e_deadtime_t,e_rawrate_t,e_rate_t,e_rawxsec_t,e_xsec_t,e_psi_t,rawrate_fit_t,rate_fit_t,rawxsec_fit_t,xsec_fit_t,e_rawrate_fit_t,e_rate_fit_t,e_rawxsec_fit_t,e_xsec_fit_t] = MakePlotArrays()
444    
445     if fit_file:
446     FitType = InputFit[print_trigger][0]
447     X0 = InputFit[print_trigger][1]
448     X1 = InputFit[print_trigger][2]
449     X2 = InputFit[print_trigger][3]
450     X3 = InputFit[print_trigger][4]
451     Chi2 = InputFit[print_trigger][5]
452 grchrist 1.14 #print str(print_trigger)+" "+str(FitType)+" "+str(X0)+" "+str(X1)+" "+str(X2)+" "+str(X3)
453     if (first_trigger):
454     print '%20s % 10s % 6s % 5s % 5s % 3s % 4s' % ('trigger', 'fit type ', 'cubic', 'quad', ' linear', ' c ', 'Chi2')
455     first_trigger=False
456     print '%20s % 10s % 2.2g % 2.2g % 2.2g % 2.2g % 2.2g' % (print_trigger, FitType, X3, X2, X1, X0, Chi2)
457     #print '{}, {}, {:02.2g}, {:02.2g}, {:02.2g}, {:02.2g} '.format(print_trigger, FitType, X0, X1, X2, X3)
458 grchrist 1.8 ## we are 2 lumis off when we start! -gets worse when we skip lumis
459     it_offset=2
460 abrinke1 1.1 for iterator in range(len(Rates[print_trigger]["rate"])):
461     if not Rates[print_trigger]["run"][iterator] in run_list:
462     continue
463     prediction = meanxsec + slopexsec * (Rates[print_trigger]["live_lumi"][iterator] - meanlumi)
464     realvalue = Rates[print_trigger]["xsec"][iterator]
465 grchrist 1.8
466 grchrist 1.11
467     if pass_cuts(data_clean, realvalue, prediction, meanxsec, Rates, print_trigger, iterator, num_ls,LumiPageInfo,SubSystemOff,max_dt,print_info, trig_list):
468 grchrist 1.8
469     if num_ls==1:
470     ##fit is 2 ls ahead of real rate
471     fit_iterator=iterator+it_offset
472     if fit_iterator>(len(Rates[print_trigger]["rate"])-1):
473     ##don't let fit_iterator go above the length of the array
474     fit_iterator=iterator
475 grchrist 1.12 else:
476     fit_iterator=iterator
477    
478 abrinke1 1.1 run_t.append(Rates[print_trigger]["run"][iterator])
479     ls_t.append(Rates[print_trigger]["ls"][iterator])
480     ps_t.append(Rates[print_trigger]["ps"][iterator])
481     inst_t.append(Rates[print_trigger]["inst_lumi"][iterator])
482 grchrist 1.8 live_t.append(Rates[print_trigger]["live_lumi"][fit_iterator])
483     delivered_t.append(Rates[print_trigger]["delivered_lumi"][fit_iterator])
484     deadtime_t.append(Rates[print_trigger]["deadtime"][fit_iterator])
485 abrinke1 1.1 rawrate_t.append(Rates[print_trigger]["rawrate"][iterator])
486     rate_t.append(Rates[print_trigger]["rate"][iterator])
487     rawxsec_t.append(Rates[print_trigger]["rawxsec"][iterator])
488     xsec_t.append(Rates[print_trigger]["xsec"][iterator])
489 abrinke1 1.5 psi_t.append(Rates[print_trigger]["psi"][iterator])
490 abrinke1 1.1
491     e_run_t.append(0.0)
492     e_ls_t.append(0.0)
493     e_ps_t.append(0.0)
494     e_inst_t.append(14.14)
495     e_live_t.append(14.14)
496     e_delivered_t.append(14.14)
497 abrinke1 1.2 e_deadtime_t.append(0.01)
498 abrinke1 1.1 e_rawrate_t.append(math.sqrt(Rates[print_trigger]["rawrate"][iterator]/(num_ls*23.3)))
499     e_rate_t.append(Rates[print_trigger]["ps"][iterator]*math.sqrt(Rates[print_trigger]["rawrate"][iterator]/(num_ls*23.3)))
500 abrinke1 1.5 e_psi_t.append(0.0)
501 abrinke1 1.2 if live_t[-1] == 0:
502     e_rawxsec_t.append(0)
503     e_xsec_t.append(0)
504     else:
505 grchrist 1.16 try:
506     e_rawxsec_t.append(math.sqrt(Rates[print_trigger]["rawrate"][iterator]/(num_ls*23.3))/Rates[print_trigger]["live_lumi"][iterator])
507     e_xsec_t.append(Rates[print_trigger]["ps"][iterator]*math.sqrt(Rates[print_trigger]["rawrate"][iterator]/(num_ls*23.3))/Rates[print_trigger]["live_lumi"][iterator])
508     except:
509     e_rawxsec_t.append(0.)
510     e_xsec_t.append(0.)
511 abrinke1 1.5 if fit_file:
512     if FitType == "expo":
513     rate_prediction = X0 + X1*math.exp(X2*delivered_t[-1])
514     else:
515     rate_prediction = X0 + X1*delivered_t[-1] + X2*delivered_t[-1]*delivered_t[-1] + X3*delivered_t[-1]*delivered_t[-1]*delivered_t[-1]
516     ## if rate_t[-1] < 0.7 * rate_prediction or rate_t[-1] > 1.4 * rate_prediction:
517     ## print str(run_t[-1])+" "+str(ls_t[-1])+" "+str(print_trigger)+" "+str(ps_t[-1])+" "+str(deadtime_t[-1])+" "+str(rate_prediction)+" "+str(rate_t[-1])+" "+str(rawrate_t[-1])
518    
519 abrinke1 1.2 if live_t[-1] == 0:
520 abrinke1 1.5 rawrate_fit_t.append(0)
521     rate_fit_t.append(0)
522 abrinke1 1.2 rawxsec_fit_t.append(0)
523     xsec_fit_t.append(0)
524 abrinke1 1.7 e_rawrate_fit_t.append(0)
525     e_rate_fit_t.append(math.sqrt(Chi2))
526     e_rawxsec_fit_t.append(0)
527     e_xsec_fit_t.append(0)
528 grchrist 1.11 #print "live_t=0", ls_t[-1], rawrate_fit_t[-1]
529 abrinke1 1.2 else:
530 grchrist 1.11 if ps_t[-1]>0.0:
531     rawrate_fit_t.append(rate_prediction*(1.0-deadtime_t[-1])/(ps_t[-1]))
532     else:
533     rawrate_fit_t.append(0.0)
534    
535 abrinke1 1.5 rate_fit_t.append(rate_prediction)
536     rawxsec_fit_t.append(rawrate_fit_t[-1]/live_t[-1])
537     xsec_fit_t.append(rate_prediction*(1.0-deadtime_t[-1])/live_t[-1])
538 abrinke1 1.7 e_rawrate_fit_t.append(math.sqrt(Chi2)*rawrate_fit_t[-1]/rate_fit_t[-1])
539     e_rate_fit_t.append(math.sqrt(Chi2))
540     e_rawxsec_fit_t.append(math.sqrt(Chi2)*rawxsec_fit_t[-1]/rate_fit_t[-1])
541     e_xsec_fit_t.append(math.sqrt(Chi2)*xsec_fit_t[-1]/rate_fit_t[-1])
542 grchrist 1.11 #print "live_t>0", ls_t[-1], rawrate_fit_t[-1]
543 abrinke1 1.5
544 grchrist 1.16 ##print iterator, fit_iterator, "ls=",ls_t[-1],"rate=",round(rawrate_t[-1],2), "deadtime=",round(deadtime_t[-1],2),"rawrate_fit=",round(rawrate_fit_t[-1],2),"max it=",len(Rates[print_trigger]["rate"])
545    
546 grchrist 1.8
547    
548 abrinke1 1.5 else: ##If the data point does not pass the data_clean filter
549 grchrist 1.11 #print "not passed", iterator, ls_t[-1], rawrate_fit_t[-1]
550 abrinke1 1.1 if debug_print:
551     print str(print_trigger)+" has xsec "+str(round(Rates[print_trigger]["xsec"][iterator],6))+" at lumi "+str(round(Rates[print_trigger]["live_lumi"][iterator],2))+" where the expected value is "+str(prediction)
552    
553 abrinke1 1.5 ## End "for iterator in range(len(Rates[print_trigger]["rate"])):" loop
554 abrinke1 1.1
555 abrinke1 1.5 AllPlotArrays = [run_t,ls_t,ps_t,inst_t,live_t,delivered_t,deadtime_t,rawrate_t,rate_t,rawxsec_t,xsec_t,psi_t,e_run_t,e_ls_t,e_ps_t,e_inst_t,e_live_t,e_delivered_t,e_deadtime_t,e_rawrate_t,e_rate_t,e_rawxsec_t,e_xsec_t,e_psi_t,rawrate_fit_t,rate_fit_t,rawxsec_fit_t,xsec_fit_t,e_rawrate_fit_t,e_rate_fit_t,e_rawxsec_fit_t,e_xsec_fit_t]
556     [VX, VXE, VY, VYE, VF, VFE] = GetVXVY(plot_properties, fit_file, AllPlotArrays)
557 abrinke1 1.1
558 abrinke1 1.5 if save_root or save_png:
559     c1 = TCanvas(str(varX),str(varY))
560     c1.SetName(str(print_trigger)+"_"+str(varY)+"_vs_"+str(varX))
561    
562     gr1 = TGraphErrors(len(VX), VX, VY, VXE, VYE)
563     gr1.SetName("Graph_"+str(print_trigger)+"_"+str(varY)+"_vs_"+str(varX))
564     gr1.GetXaxis().SetTitle(varX)
565     gr1.GetYaxis().SetTitle(varY)
566     gr1.SetTitle(str(print_trigger))
567     gr1.SetMinimum(0)
568     gr1.SetMaximum(1.2*max(VY))
569     #gr1.GetXaxis().SetLimits(min(VX)-0.2*max(VX),1.2*max(VX))
570     gr1.GetXaxis().SetLimits(0,1.2*max(VX))
571     gr1.SetMarkerStyle(8)
572     if fit_file:
573     gr1.SetMarkerSize(0.8)
574     else:
575     gr1.SetMarkerSize(0.5)
576     gr1.SetMarkerColor(2)
577    
578     if fit_file:
579     gr3 = TGraphErrors(len(VX), VX, VF, VXE, VFE)
580     gr3.SetMarkerStyle(8)
581     gr3.SetMarkerSize(0.4)
582     gr3.SetMarkerColor(4)
583     gr3.SetFillColor(4)
584     gr3.SetFillStyle(3003)
585 abrinke1 1.1
586 abrinke1 1.5 if do_fit:
587     if "rate" in varY:
588     f1a = TF1("f1a","pol2",0,8000)
589     f1a.SetLineColor(4)
590     f1a.SetLineWidth(2)
591     f1a.SetParLimits(0,0,0.2*(sum(VY)/len(VY))+0.8*min(VY))
592     f1a.SetParLimits(1,0,2.0*max(VY)/(max(VX)*max(VX)))
593     #gr1.Fit("f1a","B","Q")
594     gr1.Fit("f1a","Q","rob=0.90")
595    
596     f1b = 0
597     f1c = 0
598     if True:
599     f1b = TF1("f1b","pol3",0,8000)
600     f1b.SetLineColor(2)
601     f1b.SetLineWidth(2)
602     f1b.SetParLimits(0,0,0.2*(sum(VY)/len(VY))+0.8*min(VY))
603     f1b.SetParLimits(1,0,f1a.GetParameter(1)+0.0000001)
604     f1b.SetParLimits(2,0,f1a.GetParameter(2)+0.0000000001)
605     f1b.SetParLimits(3,0,2.0*max(VY)/(max(VX)*max(VX)*max(VX)))
606     gr1.Fit("f1b","Q","rob=0.90")
607     #if f1b.GetChisquare()/f1b.GetNDF() < f1a.GetChisquare()/f1a.GetNDF():
608 grchrist 1.15 #print "X0 = "+str(f1a.GetParameter(0))+" X1 = "+str(f1a.GetParameter(1))+" X2 = "+str(f1a.GetParameter(2))
609     #print str(print_trigger)+" f1a Chi2 = "+str(10*f1a.GetChisquare()*math.sqrt(len(VY))/(math.sqrt(sum(VY))*num_ls*f1a.GetNDF()))+", f1b Chi2 = "+str(10*f1b.GetChisquare()*math.sqrt(len(VY))/(math.sqrt(sum(VY))*num_ls*f1b.GetNDF()))
610     #print "X0 = "+str(f1b.GetParameter(0))+" X1 = "+str(f1b.GetParameter(1))+" X2 = "+str(f1b.GetParameter(2))+" X3 = "+str(f1b.GetParameter(3))
611     if (first_trigger):
612     print "len(VX)=",len(VX), "len(VY)=",len(VY)
613 grchrist 1.16 print '%-60s %4s x0 x1 x2 x3 chi2 ndf chi2/ndf' % ('trigger', 'type')
614 grchrist 1.15
615     first_trigger=False
616    
617    
618 abrinke1 1.1
619 abrinke1 1.5 f1c = TF1("f1c","[0]+[1]*expo(2)",0,8000)
620     f1c.SetLineColor(3)
621     f1c.SetLineWidth(2)
622     f1c.SetParLimits(0,0,0.2*(sum(VY)/len(VY))+0.8*min(VY))
623     f1c.SetParLimits(1,max(VY)/math.exp(10.0),max(VY)/math.exp(2.0))
624     f1c.SetParLimits(2,0.0,0.0000000001)
625     f1c.SetParLimits(3,2.0/max(VX),10.0/max(VX))
626 grchrist 1.15 #print str(max(VY)/math.exp(2.0))+" "+str(10.0/max(VX))
627 abrinke1 1.5 gr1.Fit("f1c","Q","rob=0.90")
628     #if f1c.GetChisquare()/f1c.GetNDF() < f1a.GetChisquare()/f1a.GetNDF():
629 grchrist 1.15 #print str(print_trigger)+" f1a Chi2 = "+str(10*f1a.GetChisquare()*math.sqrt(len(VY))/(math.sqrt(sum(VY))*num_ls*f1a.GetNDF()))+", f1c Chi2 = "+str(10*f1c.GetChisquare()*math.sqrt(len(VY))/(math.sqrt(sum(VY))*num_ls*f1c.GetNDF()))
630     #print "X0 = "+str(f1c.GetParameter(0))+" X1 = "+str(f1c.GetParameter(1))+" X2 = "+str(f1c.GetParameter(2))+" X3 = "+str(f1c.GetParameter(3))
631    
632    
633    
634    
635    
636     if (f1c.GetChisquare()/f1c.GetNDF() < f1b.GetChisquare()/f1b.GetNDF() and f1c.GetChisquare()/f1c.GetNDF() < f1a.GetChisquare()/f1a.GetNDF()):
637 grchrist 1.16 print '%-60s expo % .2f+/-%.2f % .2e+/-%.1e % .2e+/-%.1e % .2e+/-%.1e %7.2f %4.0f %5.3f ' % (print_trigger, f1c.GetParameter(0), f1c.GetParError(0), f1c.GetParameter(1), f1c.GetParError(1), 0 , 0 , 0 , 0 , f1c.GetChisquare(), f1c.GetNDF(), f1c.GetChisquare()/f1c.GetNDF())
638 grchrist 1.15 elif (f1b.GetChisquare()/f1b.GetNDF() < f1a.GetChisquare()/f1a.GetNDF()):
639 grchrist 1.16 print '%-60s cube % .2f+/-%.2f % .2e+/-%.1e % .2e+/-%.1e % .2e+/-%.1e %7.2f %4.0f %5.3f ' % (print_trigger, f1b.GetParameter(0), f1b.GetParError(0), f1b.GetParameter(1), f1b.GetParError(1), f1b.GetParameter(2), f1b.GetParError(2), f1b.GetParameter(3), f1b.GetParError(3), f1b.GetChisquare(), f1b.GetNDF(), f1b.GetChisquare()/f1b.GetNDF())
640 grchrist 1.15 else:
641 grchrist 1.16 print '%-60s quad % .2f+/-%.2f % .2e+/-%.1e % .2e+/-%.1e % .2e+/-%.1e %7.2f %4.0f %5.3f ' % (print_trigger, f1a.GetParameter(0), f1a.GetParError(0), f1a.GetParameter(1), f1a.GetParError(1), f1a.GetParameter(2), f1a.GetParError(2), 0 , 0 , f1a.GetChisquare(), f1a.GetNDF(), f1a.GetChisquare()/f1a.GetNDF())
642 abrinke1 1.2
643 grchrist 1.15
644    
645 abrinke1 1.5 else: ##If this is not a rate plot
646     f1a = TF1("f1a","pol1",0,8000)
647     f1a.SetLineColor(4)
648     f1a.SetLineWidth(2)
649     if "xsec" in varY:
650     f1a.SetParLimits(0,0,meanxsec*1.5)
651     if slopexsec > 0:
652     f1a.SetParLimits(1,0,max(VY)/max(VX))
653     else:
654     f1a.SetParLimits(1,2*slopexsec,-2*slopexsec)
655 abrinke1 1.1 else:
656 abrinke1 1.5 f1a.SetParLimits(0,-1000,1000)
657     gr1.Fit("f1a","Q","rob=0.80")
658 abrinke1 1.1
659 abrinke1 1.5 if save_root or save_png:
660     gr1.Draw("APZ")
661     ## ##Option to draw stats box
662     ## p1 = TPaveStats()
663     ## p1 = gr1.GetListOfFunctions().FindObject("stats")
664     ## print p1
665     ## gr1.PaintStats(f1b).Draw("same")
666     if fit_file:
667     gr3.Draw("P3")
668     if do_fit:
669     f1a.Draw("same")
670     try:
671     f1b.Draw("same")
672     f1c.Draw("same")
673     except:
674     True
675     c1.Update()
676     if save_root:
677     myfile = TFile( RootFile, 'UPDATE' )
678     c1.Write()
679     myfile.Close()
680     if save_png:
681     c1.SaveAs(str(print_trigger)+"_"+str(varY)+"_vs_"+str(varX)+".png")
682 abrinke1 1.1
683    
684     if print_table or save_fits:
685 abrinke1 1.5 if not do_fit:
686     print "Can't have save_fits = True and do_fit = False"
687     continue
688     if f1c.GetChisquare()/f1c.GetNDF() < 0.95*f1a.GetChisquare()/f1a.GetNDF() and f1c.GetChisquare()/f1c.GetNDF() < 0.95*f1b.GetChisquare()/f1b.GetNDF():
689 grchrist 1.16 OutputFit[print_trigger] = ["expo", f1c.GetParameter(0), f1c.GetParameter(1), f1c.GetParameter(3), 0.0, f1c.GetChisquare()/f1c.GetNDF(), meanrawrate, f1c.GetParError(0), f1c.GetParError(1), f1c.GetParError(2), f1c.GetParError(3)]
690 abrinke1 1.5 elif f1b.GetChisquare()/f1b.GetNDF() < 0.95*f1a.GetChisquare()/f1a.GetNDF():
691 grchrist 1.16 OutputFit[print_trigger] = ["poly", f1b.GetParameter(0), f1b.GetParameter(1), f1b.GetParameter(2), f1b.GetParameter(3), f1b.GetChisquare()/f1b.GetNDF(), meanrawrate,f1b.GetParError(0), f1b.GetParError(1), f1b.GetParError(2), f1b.GetParError(3)]
692 abrinke1 1.5 else:
693 grchrist 1.16 OutputFit[print_trigger] = ["poly", f1a.GetParameter(0), f1a.GetParameter(1), f1a.GetParameter(2), 0.0, f1a.GetChisquare()/f1a.GetNDF(), meanrawrate, f1a.GetParError(0), f1a.GetParError(1), f1a.GetParError(2), 0.0]
694 abrinke1 1.1
695     if save_root:
696     print "Output root file is "+str(RootFile)
697    
698     if save_fits:
699     FitNameTemplate = "Fits/2011/Fit_%s_%sLS_Run%sto%s.pkl"
700     FitFile = FitNameTemplate % (trig_name, num_ls, min_run, max_run)
701     if os.path.exists(FitFile):
702     os.remove(FitFile)
703 abrinke1 1.5 FitOutputFile = open(FitFile, 'wb')
704     pickle.dump(OutputFit, FitOutputFile, 2)
705     FitOutputFile.close()
706 abrinke1 1.7 print "Output fit file is "+str(FitFile)
707 abrinke1 1.1
708     if print_table:
709 abrinke1 1.5 print "The expo fit is of the form p0+p1*e^(p2*x), poly is p0+(p1/10^3)*x+(p2/10^6)*x^2+(p3/10^9)*x^3, where x is Deliv. Lumi."
710     print '%60s%10s%10s%10s%10s%10s%10s%10s' % ("Trig", "fit", "p0", "p1", "p2", "p3", "Chi2", "Av raw")
711     for print_trigger in OutputFit:
712     _trigger = (print_trigger[:56] + '...') if len(print_trigger) > 59 else print_trigger
713 abrinke1 1.1 try:
714 abrinke1 1.5 if OutputFit[print_trigger][0] == "poly":
715     print '%60s%10s%10s%10s%10s%10s%10s' % (_trigger, OutputFit[print_trigger][0], round(OutputFit[print_trigger][1],3), round(OutputFit[print_trigger][2],6)*1000, round(OutputFit[print_trigger][3],9)*1000000, round(OutputFit[print_trigger][4],12)*1000000000, round(OutputFit[print_trigger][5],2), round(OutputFit[print_trigger][6],3))
716     else:
717     print '%60s%10s%10s%10s%10s%10s%10s' % (_trigger, OutputFit[print_trigger][0], OutputFit[print_trigger][1], OutputFit[print_trigger][2], OutputFit[print_trigger][3], OutputFit[print_trigger][4], round(OutputFit[print_trigger][5],2), round(OutputFit[print_trigger][6],3))
718 abrinke1 1.1 except:
719     print str(print_trigger)+" is somehow broken"
720    
721 abrinke1 1.5
722     ############# SUPPORTING FUNCTIONS ################
723    
724    
725     def GetJSON(json_file):
726    
727     input_file = open(json_file)
728     file_content = input_file.read()
729     inputRange = selectionParser(file_content)
730     JSON = inputRange.runsandls()
731     return JSON
732     ##JSON is an array: JSON[run_number] = [1st ls, 2nd ls, 3rd ls ... nth ls]
733    
734     def MakePlotArrays():
735     run_t = array.array('f')
736     ls_t = array.array('f')
737     ps_t = array.array('f')
738     inst_t = array.array('f')
739     live_t = array.array('f')
740     delivered_t = array.array('f')
741     deadtime_t = array.array('f')
742     rawrate_t = array.array('f')
743     rate_t = array.array('f')
744     rawxsec_t = array.array('f')
745     xsec_t = array.array('f')
746     psi_t = array.array('f')
747    
748     e_run_t = array.array('f')
749     e_ls_t = array.array('f')
750     e_ps_t = array.array('f')
751     e_inst_t = array.array('f')
752     e_live_t = array.array('f')
753     e_delivered_t = array.array('f')
754     e_deadtime_t = array.array('f')
755     e_rawrate_t = array.array('f')
756     e_rate_t = array.array('f')
757     e_rawxsec_t = array.array('f')
758     e_xsec_t = array.array('f')
759     e_psi_t = array.array('f')
760    
761     rawrate_fit_t = array.array('f')
762     rate_fit_t = array.array('f')
763     rawxsec_fit_t = array.array('f')
764     xsec_fit_t = array.array('f')
765     e_rawrate_fit_t = array.array('f')
766     e_rate_fit_t = array.array('f')
767     e_rawxsec_fit_t = array.array('f')
768     e_xsec_fit_t = array.array('f')
769    
770     return [run_t,ls_t,ps_t,inst_t,live_t,delivered_t,deadtime_t,rawrate_t,rate_t,rawxsec_t,xsec_t,psi_t,e_run_t,e_ls_t,e_ps_t,e_inst_t,e_live_t,e_delivered_t,e_deadtime_t,e_rawrate_t,e_rate_t,e_rawxsec_t,e_xsec_t,e_psi_t,rawrate_fit_t,rate_fit_t,rawxsec_fit_t,xsec_fit_t,e_rawrate_fit_t,e_rate_fit_t,e_rawxsec_fit_t,e_xsec_fit_t]
771    
772    
773     def GetVXVY(plot_properties, fit_file, AllPlotArrays):
774    
775     VF = "0"
776     VFE = "0"
777    
778     [run_t,ls_t,ps_t,inst_t,live_t,delivered_t,deadtime_t,rawrate_t,rate_t,rawxsec_t,xsec_t,psi_t,e_run_t,e_ls_t,e_ps_t,e_inst_t,e_live_t,e_delivered_t,e_deadtime_t,e_rawrate_t,e_rate_t,e_rawxsec_t,e_xsec_t,e_psi_t,rawrate_fit_t,rate_fit_t,rawxsec_fit_t,xsec_fit_t,e_rawrate_fit_t,e_rate_fit_t,e_rawxsec_fit_t,e_xsec_fit_t] = AllPlotArrays
779     for varX, varY, do_fit, save_root, save_png, fit_file in plot_properties:
780     if varX == "run":
781     VX = run_t
782     VXE = run_t_e
783     elif varX == "ls":
784     VX = ls_t
785     VXE = e_ls_t
786     elif varX == "ps":
787     VX = ps_t
788     VXE = e_ps_t
789     elif varX == "inst":
790     VX = inst_t
791     VXE = e_inst_t
792     elif varX == "live":
793     VX = live_t
794     VXE = e_live_t
795     elif varX == "delivered":
796     VX = delivered_t
797     VXE = e_delivered_t
798     elif varX == "deadtime":
799     VX = deadtime_t
800     VXE = e_deadtime_t
801     elif varX == "rawrate":
802     VX = rawrate_t
803     VXE = e_rawrate_t
804     elif varX == "rate":
805     VX = rate_t
806     VXE = e_rate_t
807     elif varX == "rawxsec":
808     VX = rawxsec_t
809     VXE = e_rawxsec_t
810     elif varX == "xsec":
811     VX = xsec_t
812     VXE = e_xsec_t
813     elif varX == "psi":
814     VX = psi_t
815     VXE = e_psi_t
816     else:
817     print "No valid variable entered for X"
818     continue
819     if varY == "run":
820     VY = run_t
821     VYE = run_t_e
822     elif varY == "ls":
823     VY = ls_t
824     VYE = e_ls_t
825     elif varY == "ps":
826     VY = ps_t
827     VYE = e_ps_t
828     elif varY == "inst":
829     VY = inst_t
830     VYE = e_inst_t
831     elif varY == "live":
832     VY = live_t
833     VYE = e_live_t
834     elif varY == "delivered":
835     VY = delivered_t
836     VYE = e_delivered_t
837     elif varY == "deadtime":
838     VY = deadtime_t
839     VYE = e_deadtime_t
840     elif varY == "rawrate":
841     VY = rawrate_t
842     VYE = e_rawrate_t
843     if fit_file:
844     VF = rawrate_fit_t
845     VFE = e_rawrate_fit_t
846     elif varY == "rate":
847     VY = rate_t
848     VYE = e_rate_t
849     if fit_file:
850     VF = rate_fit_t
851     VFE = e_rate_fit_t
852     elif varY == "rawxsec":
853     VY = rawxsec_t
854     VYE = e_rawxsec_t
855     if fit_file:
856     VF = rawxsec_fit_t
857     VFE = e_rawxsec_fit_t
858     elif varY == "xsec":
859     VY = xsec_t
860     VYE = e_xsec_t
861     if fit_file:
862     VF = xsec_fit_t
863     VFE = e_xsec_fit_t
864     elif varY == "psi":
865     VY = psi_t
866     VYE = e_psi_t
867     else:
868     print "No valid variable entered for Y"
869     continue
870    
871     return [VX, VXE, VY, VYE, VF, VFE]
872    
873 grchrist 1.11 def pass_cuts(data_clean, realvalue, prediction, meanxsec, Rates, print_trigger, iterator, num_ls,LumiPageInfo,SubSystemOff, max_dt, print_info, trig_list):
874 grchrist 1.8 it_offset=2
875 grchrist 1.15 Passed=True
876     subsystemfailed=[]
877 grchrist 1.11
878 grchrist 1.8 if num_ls==1:
879     ##fit is 2 ls ahead of real rate
880     fit_iterator=iterator+it_offset
881     if fit_iterator>(len(Rates[print_trigger]["rate"])-1):
882     ##don't let fit_iterator go above the length of the array
883     fit_iterator=iterator
884    
885 grchrist 1.10 LS=Rates[print_trigger]["ls"][iterator]
886     #print "ls=",LS,
887     LSRange=LumiPageInfo[LS]["LSRange"]
888     #print LSRange,
889     LS2=LSRange[-1]
890     #LS2=LSRange.pop()
891     #print "LS2=",LS2
892    
893     #print LumiPageInfo[LS]
894     lumidict={}
895     lumidict=LumiPageInfo[LS]
896 grchrist 1.15
897 grchrist 1.11
898    
899    
900    
901     if print_info:
902     if (iterator==0 and print_trigger==trig_list[0]):
903     print '%10s%10s%10s%10s%10s%10s%10s%15s%20s' % ("Status", "Run", "LS", "Physics", "Active", "Deadtime", " MaxDeadTime", " Passed all subsystems?", " List of Subsystems failed")
904    
905     ## if SubSystemOff["All"]:
906     ## for keys in LumiPageInfo[LS]:
907     ## #print LS, keys, LumiPageInfo[LS][keys]
908     ## if not LumiPageInfo[LS][keys]:
909     ## Passed=False
910     ## subsystemfailed.append(keys)
911     ## break
912     ## else:
913     if SubSystemOff["Mu"] or SubSystemOff["All"]:
914     if not (LumiPageInfo[LS]["rpc"] or LumiPageInfo[LS]["dt0"] or LumiPageInfo[LS]["dtp"] or LumiPageInfo[LS]["dtm"] or LumiPageInfo["cscp"] or LumiPageInfo["cscm"]):
915     Passed=False
916     subsystemfailed.append("Mu")
917     if SubSystemOff["HCal"] or SubSystemOff["All"]:
918     if not (LumiPageInfo[LS]["hbhea"] and LumiPageInfo[LS]["hbheb"] and LumiPageInfo[LS]["hbhec"]):
919     Passed=False
920     subsystemfailed.append("HCal")
921     if (SubSystemOff["EndCap"] or SubSystemOff["All"]) and not (LumiPageInfo[LS]["hf"]):
922     Passed=False
923     subsystemfailed.append("HCal-EndCap")
924     if SubSystemOff["ECal"] or SubSystemOff["All"]:
925     if not (LumiPageInfo[LS]["ebp"] and LumiPageInfo[LS]["ebm"]):
926     Passed=False
927     subsystemfailed.append("ECal")
928     if (SubSystemOff["EndCap"] or SubSystemOff["All"]) and not (LumiPageInfo[LS]["eep"] and LumiPageInfo[LS]["eem"] and LumiPageInfo[LS]["esp"] or LumiPageInfo[LS]["esm"]):
929     Passed=False
930     subsystemfailed.append("ECal-EndCap")
931     if SubSystemOff["Tracker"] or SubSystemOff["All"]:
932     if not (LumiPageInfo[LS]["tob"] and LumiPageInfo[LS]["tibtid"] and LumiPageInfo[LS]["bpix"] and LumiPageInfo[LS]["fpix"]):
933     Passed=False
934     subsystemfailed.append("Tracker")
935     if (SubSystemOff["EndCap"] or SubSystemOff["All"]) and not (LumiPageInfo[LS]["tecp"] and LumiPageInfo[LS]["tecm"]):
936     Passed=False
937     subsystemfailed.append("Tracker-EndCap")
938     if SubSystemOff["Beam"] or SubSystemOff["All"]:
939     if not(LumiPageInfo[LS]["b1pres"] and LumiPageInfo[LS]["b2pres"] and LumiPageInfo[LS]["b1stab"] and LumiPageInfo[LS]["b2stab"]):
940     Passed=False
941     subsystemfailed.append("Beam")
942 grchrist 1.12 else:
943     fit_iterator=iterator
944     Passed=True
945 grchrist 1.10
946 grchrist 1.11 #print "LS",LS, Passed, round(Rates[print_trigger]["deadtime"][fit_iterator],2), max_dt
947 grchrist 1.10
948 grchrist 1.8 if not data_clean or (
949 grchrist 1.11 ## (
950 grchrist 1.8
951 grchrist 1.11 ## (
952     ## realvalue > 0.4*prediction
953     ## and realvalue < 2.5*prediction
954     ## )
955 grchrist 1.8
956 grchrist 1.11 ## or
957 grchrist 1.8
958 grchrist 1.11 ## (
959     ## realvalue > 0.4*meanxsec
960     ## and realvalue < 2.5*meanxsec
961     ## )
962 grchrist 1.8
963 grchrist 1.11 ## or prediction < 0
964     ## )
965 grchrist 1.8
966 grchrist 1.11 Rates[print_trigger]["physics"][iterator] == 1
967 grchrist 1.8 and Rates[print_trigger]["active"][iterator] == 1
968 grchrist 1.11 and Rates[print_trigger]["deadtime"][fit_iterator] < max_dt
969     #and Rates[print_trigger]["psi"][iterator] > 0
970 grchrist 1.10 and Passed
971 grchrist 1.8 ):
972 grchrist 1.11 #print LS, "True"
973 grchrist 1.8 return True
974     else:
975 grchrist 1.11
976 grchrist 1.13 if (print_info and print_trigger==trig_list[0] and num_ls==1):
977 grchrist 1.11
978     print '%10s%10s%10s%10s%10s%10s%10s%15s%20s' % ("Failed", Rates[print_trigger]["run"][iterator], LS, Rates[print_trigger]["physics"][iterator], Rates[print_trigger]["active"][iterator], round(Rates[print_trigger]["deadtime"][fit_iterator],2), max_dt, Passed, subsystemfailed)
979 grchrist 1.15 ##elif(print_info and print_trigger==trig_list[0]):
980     ## print '%10s%10s%10s%10s%10s%10s%10s%15s%20s' % ("Failed", Rates[print_trigger]["run"][iterator], LumiPageInfo[Rates[print_trigger]["ls"][iterator]]["LSRange"], Rates[print_trigger]["physics"][iterator], Rates[print_trigger]["active"][iterator], round(Rates[print_trigger]["deadtime"][fit_iterator],2), max_dt, Passed, subsystemfailed)
981 grchrist 1.8 return False
982 abrinke1 1.5
983 grchrist 1.10
984     #### LumiRangeGreens ####
985     ####inputs: RefMoreLumiArray --dict with lumi page info in LS by LS blocks,
986     #### LRange --list range over lumis,
987     #### nls --number of lumisections
988     #### RefRunNum --run number
989     ####
990     ####outputs RangeMoreLumi --lumi page info in dict LSRange blocks with lumi, added items Run and LSRange
991     def LumiRangeGreens(RefMoreLumiArray,LSRange,nls,RefRunNum):
992 grchrist 1.9
993     RangeMoreLumi={}
994     for keys,values in RefMoreLumiArray.iteritems():
995     RangeMoreLumi[keys]=1
996 grchrist 1.10
997 grchrist 1.9 for iterator in LSRange[nls]:
998     for keys, values in RefMoreLumiArray.iteritems():
999     if RefMoreLumiArray[keys][iterator]==0:
1000     RangeMoreLumi[keys]=0
1001 grchrist 1.10 RangeMoreLumi['LSRange']=LSRange[nls]
1002     RangeMoreLumi['Run']=RefRunNum
1003    
1004 grchrist 1.9 return RangeMoreLumi
1005    
1006 grchrist 1.10 #### CheckLumis ####
1007     ####inputs:
1008     #### PageLumiInfo --dict of LS with dict of some lumipage info
1009     #### Rates --dict of triggernames with dict of info
1010     def checkLS(Rates, PageLumiInfo,trig_list):
1011     rateslumis=Rates[trig_list[-1]]["ls"]
1012     keys=PageLumiInfo.keys()
1013     print "lumi run=",PageLumiInfo[keys[-1]]["Run"]
1014     ll=0
1015     for ls in keys:
1016     print ls,rateslumis[ll]
1017     ll=ll+1
1018     return False
1019    
1020    
1021 grchrist 1.9
1022 abrinke1 1.1 if __name__=='__main__':
1023     main()