ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/RateMonShiftTool_dev/DatabaseRatePredictor.py
Revision: 1.19
Committed: Wed Mar 21 15:21:47 2012 UTC (13 years, 1 month ago) by amott
Content type: text/x-python
Branch: MAIN
CVS Tags: V00-00-16
Changes since 1.18: +0 -3 lines
Log Message:
bug fix

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