| 1 |
#!/usr/bin/env python
|
| 2 |
|
| 3 |
from DatabaseParser import *
|
| 4 |
from GetListOfRuns import *
|
| 5 |
import sys
|
| 6 |
import os
|
| 7 |
from numpy import *
|
| 8 |
import pickle
|
| 9 |
import getopt
|
| 10 |
from StreamMonitor import StreamMonitor
|
| 11 |
from itertools import groupby
|
| 12 |
from operator import itemgetter
|
| 13 |
from collections import deque
|
| 14 |
|
| 15 |
from ROOT import gROOT, TCanvas, TF1, TGraph, TGraphErrors, TPaveStats, gPad, gStyle
|
| 16 |
from ROOT import TFile, TPaveText, TBrowser
|
| 17 |
from ROOT import gBenchmark
|
| 18 |
import array
|
| 19 |
import math
|
| 20 |
from ReadConfig import RateMonConfig
|
| 21 |
from TablePrint import *
|
| 22 |
from selectionParser import selectionParser
|
| 23 |
|
| 24 |
def usage():
|
| 25 |
print sys.argv[0]+" [options] <list of runs>"
|
| 26 |
print "This script is used to generate fits and do secondary shifter validation"
|
| 27 |
print "For more information, see https://twiki.cern.ch/twiki/bin/view/CMS/RateMonitoringScriptWithReferenceComparison"
|
| 28 |
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"
|
| 29 |
print " be careful with using ranges (c-d), it is highly recommended to use a JSON in this case"
|
| 30 |
print "options: "
|
| 31 |
print "--makeFits run in fit making mode"
|
| 32 |
print "--secondary run in secondary shifter mode"
|
| 33 |
print "--fitFile=<path> path to the fit file"
|
| 34 |
print "--json=<path> path to the JSON file"
|
| 35 |
print "--TriggerList=<path> path to the trigger list (without versions!)"
|
| 36 |
print "--AllTriggers Run for all triggers instead of specifying a trigger list"
|
| 37 |
print "--maxdt=<max deadtime> Mask LS above max deadtime threshold"
|
| 38 |
print "--All Mask LS with any red LS on WBM LS page (not inc castor zdc etc)"
|
| 39 |
print "--Mu Mask LS with Mu (RPC, DT+, DT-, DT0, CSC+ and CSC-) off"
|
| 40 |
print "--HCal Mask LS with HCal barrel off"
|
| 41 |
print "--Tracker Mask LS with Tracker barrel off"
|
| 42 |
print "--ECal Mask LS with ECal barrel off"
|
| 43 |
print "--EndCap Mask LS with EndCap sys off, used in combination with other subsys"
|
| 44 |
print "--Beam Mask LS with Beam off"
|
| 45 |
print "--UseVersionNumbers Don't ignore path version numbers"
|
| 46 |
print "--linear Force linear fits"
|
| 47 |
print "--inst Make fits using instantaneous luminosity instead of delivered"
|
| 48 |
print "--write Writes fit info into csv, for ranking nonlinear triggers"
|
| 49 |
|
| 50 |
class Modes:
|
| 51 |
none,fits,secondary = range(3)
|
| 52 |
|
| 53 |
def pickYear():
|
| 54 |
global thisyear
|
| 55 |
thisyear="2012"
|
| 56 |
print "Year set to ",thisyear
|
| 57 |
|
| 58 |
def main():
|
| 59 |
gROOT.SetBatch(True)
|
| 60 |
try:
|
| 61 |
##set year to 2012
|
| 62 |
pickYear()
|
| 63 |
|
| 64 |
try:
|
| 65 |
opt, args = getopt.getopt(sys.argv[1:],"",["makeFits","secondary","fitFile=","json=","TriggerList=","maxdt=","All","Mu","HCal","Tracker","ECal","EndCap","Beam","UseVersionNumbers","linear","inst","write","AllTriggers","UsePSCol="])
|
| 66 |
|
| 67 |
except getopt.GetoptError, err:
|
| 68 |
print str(err)
|
| 69 |
usage()
|
| 70 |
sys.exit(2)
|
| 71 |
|
| 72 |
##### RUN LIST ########
|
| 73 |
run_list=[]
|
| 74 |
|
| 75 |
if len(args)<1:
|
| 76 |
inputrunlist=[]
|
| 77 |
print "No runs specified"
|
| 78 |
runinput=raw_input("Enter run range in form <run1> <run2> <run3> or <run1>-<run2>:")
|
| 79 |
inputrunlist.append(runinput)
|
| 80 |
|
| 81 |
if runinput.find(' ')!=-1:
|
| 82 |
args=runinput.split(' ')
|
| 83 |
else:
|
| 84 |
args.append(runinput)
|
| 85 |
|
| 86 |
for r in args:
|
| 87 |
if r.find('-')!=-1: # r is a run range
|
| 88 |
rrange = r.split('-')
|
| 89 |
if len(rrange)!=2:
|
| 90 |
print "Invalid run range %s" % (r,)
|
| 91 |
sys.exit(0)
|
| 92 |
try:
|
| 93 |
for rr in range(int(rrange[0]),int(rrange[1])+1):
|
| 94 |
run_list.append(rr)
|
| 95 |
except:
|
| 96 |
print "Invalid run range %s" % (r,)
|
| 97 |
sys.exit(0)
|
| 98 |
else: # r is not a run range
|
| 99 |
try:
|
| 100 |
run_list.append(int(r))
|
| 101 |
except:
|
| 102 |
print "Invalid run %s" % (r,)
|
| 103 |
|
| 104 |
##### READ CMD LINE ARGS #########
|
| 105 |
mode = Modes.none
|
| 106 |
fitFile = ""
|
| 107 |
jsonfile = ""
|
| 108 |
trig_list = []
|
| 109 |
max_dt=-1.0
|
| 110 |
subsys=-1.0
|
| 111 |
NoVersion=True
|
| 112 |
linear=False
|
| 113 |
do_inst=False
|
| 114 |
wp_bool=False
|
| 115 |
all_triggers=False
|
| 116 |
DoL1=True
|
| 117 |
UsePSCol=-1
|
| 118 |
SubSystemOff={'All':False,'Mu':False,'HCal':False,'ECal':False,'Tracker':False,'EndCap':False,'Beam':False}
|
| 119 |
for o,a in opt:
|
| 120 |
if o == "--makeFits":
|
| 121 |
mode = Modes.fits
|
| 122 |
elif o == "--secondary":
|
| 123 |
mode = Modes.secondary
|
| 124 |
elif o == "--fitFile":
|
| 125 |
fitFile = str(a)
|
| 126 |
elif o == "--json":
|
| 127 |
jsonfile = a
|
| 128 |
elif o=="--maxdt":
|
| 129 |
max_dt = float(a)
|
| 130 |
elif o=="--All":
|
| 131 |
subsys=1
|
| 132 |
SubSystemOff["All"]=True
|
| 133 |
elif o=="--Mu":
|
| 134 |
subsys=1
|
| 135 |
SubSystemOff["Mu"]=True
|
| 136 |
elif o=="--HCal":
|
| 137 |
SubSystemOff["HCal"]=True
|
| 138 |
subsys=1
|
| 139 |
elif o=="--Tracker":
|
| 140 |
SubSystemOff["Tracker"]=True
|
| 141 |
subsys=1
|
| 142 |
elif o=="--ECal":
|
| 143 |
SubSystemOff["ECal"]=True
|
| 144 |
subsys=1
|
| 145 |
elif o=="--EndCap":
|
| 146 |
SubSystemOff["EndCap"]=True
|
| 147 |
subsys=1
|
| 148 |
elif o=="--Beam":
|
| 149 |
SubSystemOff["Beam"]=True
|
| 150 |
subsys=1
|
| 151 |
elif o=="--UseVersionNumbers":
|
| 152 |
NoVersion=False
|
| 153 |
elif o=="--linear":
|
| 154 |
linear=True
|
| 155 |
elif o=="--inst":
|
| 156 |
do_inst=True
|
| 157 |
elif o=="--write":
|
| 158 |
wp_bool=True
|
| 159 |
elif o=="--AllTriggers":
|
| 160 |
all_triggers=True
|
| 161 |
elif o=="--UsePSCol":
|
| 162 |
UsePSCol=int(a)
|
| 163 |
elif o == "--TriggerList":
|
| 164 |
try:
|
| 165 |
f = open(a)
|
| 166 |
for entry in f:
|
| 167 |
if entry.startswith('#'):
|
| 168 |
continue
|
| 169 |
if entry.find(':')!=-1:
|
| 170 |
entry = entry[:entry.find(':')] ## We can point this to the existing monitor list, just remove everything after ':'!
|
| 171 |
if entry.find('#')!=-1:
|
| 172 |
entry = entry[:entry.find('#')] ## We can point this to the existing monitor list, just remove everything after ':'!
|
| 173 |
trig_list.append( entry.rstrip('\n'))
|
| 174 |
except:
|
| 175 |
print "\nInvalid Trigger List\n"
|
| 176 |
sys.exit(0)
|
| 177 |
else:
|
| 178 |
print "\nInvalid Option %s\n" % (str(o),)
|
| 179 |
usage()
|
| 180 |
sys.exit(2)
|
| 181 |
|
| 182 |
print "\n\n"
|
| 183 |
|
| 184 |
###### MODES #########
|
| 185 |
if mode == Modes.none: ## no mode specified
|
| 186 |
print "\nNo operation mode specified!\n"
|
| 187 |
modeinput=raw_input("Enter mode, --makeFits or --secondary:")
|
| 188 |
print "modeinput=",modeinput
|
| 189 |
if not (modeinput=="--makeFits" or modeinput=="--secondary"):
|
| 190 |
print "not either"
|
| 191 |
usage()
|
| 192 |
sys.exit(0)
|
| 193 |
elif modeinput == "--makeFits":
|
| 194 |
mode=Modes.fits
|
| 195 |
elif modeinput == "--secondary":
|
| 196 |
mode=Modes.secondary
|
| 197 |
else:
|
| 198 |
print "FATAL ERROR: No Mode specified"
|
| 199 |
sys.exit(0)
|
| 200 |
|
| 201 |
if mode == Modes.fits:
|
| 202 |
print "Running in Fit Making mode\n\n"
|
| 203 |
elif mode == Modes.secondary:
|
| 204 |
print "Running in Secondary Shifter mode\n\n"
|
| 205 |
else: ## should never get here, but exit if we do
|
| 206 |
print "FATAL ERROR: No Mode specified"
|
| 207 |
sys.exit(0)
|
| 208 |
|
| 209 |
if fitFile=="" and not mode==Modes.fits:
|
| 210 |
print "\nPlease specify fit file. These are available:\n"
|
| 211 |
path="Fits/%s/" % (thisyear) # insert the path to the directory of interest
|
| 212 |
dirList=os.listdir(path)
|
| 213 |
for fname in dirList:
|
| 214 |
print fname
|
| 215 |
fitFile = path+raw_input("Enter fit file in format Fit_HLT_10LS_Run176023to180252.pkl: ")
|
| 216 |
|
| 217 |
elif fitFile=="":
|
| 218 |
NoVstr=""
|
| 219 |
if NoVersion:
|
| 220 |
NoVstr="NoV_"
|
| 221 |
if not do_inst:
|
| 222 |
fitFile="Fits/%s/Fit_HLT_%s10LS_Run%sto%s.pkl" % (thisyear,NoVstr,min(run_list),max(run_list))
|
| 223 |
else:
|
| 224 |
fitFile="Fits/%s/Fit_inst_HLT_%s10LS_Run%sto%s.pkl" % (thisyear,NoVstr,min(run_list),max(run_list))
|
| 225 |
|
| 226 |
if "NoV" in fitFile:
|
| 227 |
NoVersion=True
|
| 228 |
|
| 229 |
###### TRIGGER LIST #######
|
| 230 |
if trig_list == [] and not all_triggers:
|
| 231 |
print "\nPlease specify list of triggers\n"
|
| 232 |
print "Available lists are:"
|
| 233 |
dirList=os.listdir(".")
|
| 234 |
for fname in dirList:
|
| 235 |
entry=fname
|
| 236 |
if entry.find('.')!=-1:
|
| 237 |
extension = entry[entry.find('.'):] ## We can point this to the existing monitor list, just remove everything after ':'!
|
| 238 |
if extension==".list":
|
| 239 |
print fname
|
| 240 |
trig_input=raw_input("\nEnter triggers in format HLT_IsoMu30_eta2p1 or a .list file, or enter AllTriggers to run over all triggers in the menu: ")
|
| 241 |
|
| 242 |
if trig_input.find('AllTriggers') != -1:
|
| 243 |
all_triggers = True
|
| 244 |
elif trig_input.find('.') != -1:
|
| 245 |
extension = trig_input[trig_input.find('.'):]
|
| 246 |
if extension==".list":
|
| 247 |
try:
|
| 248 |
fl=open(trig_input)
|
| 249 |
except:
|
| 250 |
print "Cannot open file"
|
| 251 |
usage()
|
| 252 |
sys.exit(0)
|
| 253 |
|
| 254 |
for line in fl:
|
| 255 |
if line.startswith('#'):
|
| 256 |
continue
|
| 257 |
if len(line)<1:
|
| 258 |
continue
|
| 259 |
if len(line)>=2:
|
| 260 |
arg=line.rstrip('\n').rstrip(' ').lstrip(' ')
|
| 261 |
trig_list.append(arg)
|
| 262 |
else:
|
| 263 |
arg=''
|
| 264 |
else:
|
| 265 |
trig_list.append(trig_input)
|
| 266 |
|
| 267 |
if jsonfile=="":
|
| 268 |
JSON=[]
|
| 269 |
else:
|
| 270 |
print "Using JSON: %s" % (jsonfile,)
|
| 271 |
JSON = GetJSON(jsonfile) ##Returns array JSON[runs][ls_list]
|
| 272 |
|
| 273 |
###### TO CREATE FITS #########
|
| 274 |
if mode == Modes.fits:
|
| 275 |
trig_name = "HLT"
|
| 276 |
num_ls = 10
|
| 277 |
physics_active_psi = True ##Requires that physics and active be on, and that the prescale column is not 0
|
| 278 |
debug_print = False
|
| 279 |
no_versions=False
|
| 280 |
min_rate = 0.0
|
| 281 |
print_table = False
|
| 282 |
data_clean = True ##Gets rid of anomalous rate points, reqires physics_active_psi (PAP) and deadtime < 20%
|
| 283 |
##plot_properties = [varX, varY, do_fit, save_root, save_png, fit_file]
|
| 284 |
if not do_inst:
|
| 285 |
plot_properties = [["delivered", "rate", True, True, False, fitFile]]
|
| 286 |
# plot_properties = [["delivered", "rawrate", True, True, False, fitFile]]
|
| 287 |
else:
|
| 288 |
plot_properties = [["inst", "rate", True, True, False, fitFile]]
|
| 289 |
|
| 290 |
masked_triggers = ["AlCa_", "DST_", "HLT_L1", "HLT_Zero", "HLT_BeamGas", "HLT_Activity", "L1_BeamGas", "L1_ZeroBias"]
|
| 291 |
save_fits = True
|
| 292 |
if max_dt==-1.0:
|
| 293 |
max_dt=0.08 ## no deadtime cutuse 2.0
|
| 294 |
force_new=True
|
| 295 |
print_info=True
|
| 296 |
if subsys==-1.0:
|
| 297 |
SubSystemOff={'All':False,'Mu':False,'HCal':False,'ECal':False,'Tracker':False,'EndCap':False,'Beam':True}
|
| 298 |
|
| 299 |
###### TO SEE RATE VS PREDICTION ########
|
| 300 |
if mode == Modes.secondary:
|
| 301 |
trig_name = "HLT"
|
| 302 |
num_ls = 1
|
| 303 |
physics_active_psi = True
|
| 304 |
debug_print = False
|
| 305 |
no_versions=False
|
| 306 |
min_rate = 0.0
|
| 307 |
print_table = False
|
| 308 |
data_clean = True
|
| 309 |
##plot_properties = [varX, varY, do_fit, save_root, save_png, fit_file]
|
| 310 |
plot_properties = [["ls", "rawrate", False, True, False,fitFile]]
|
| 311 |
## rate is calculated as: (measured rate, deadtime corrected) * prescale [prediction not dt corrected]
|
| 312 |
## rawrate is calculated as: measured rate [prediction is dt corrected]
|
| 313 |
|
| 314 |
masked_triggers = ["AlCa_", "DST_", "HLT_L1", "HLT_Zero", "HLT_BeamGas", "HLT_Activity", "L1_BeamGas", "L1_ZeroBias"]
|
| 315 |
save_fits = False
|
| 316 |
if max_dt==-1.0:
|
| 317 |
max_dt=2.0 ## no deadtime cut=2.0
|
| 318 |
force_new=True
|
| 319 |
print_info=True
|
| 320 |
if subsys==-1.0:
|
| 321 |
SubSystemOff={'All':True,'Mu':False,'HCal':False,'ECal':False,'Tracker':False,'EndCap':False,'Beam':True}
|
| 322 |
|
| 323 |
for k in SubSystemOff.iterkeys():
|
| 324 |
print k,"=",SubSystemOff[k]," ",
|
| 325 |
print " "
|
| 326 |
L1SeedChangeFit=True
|
| 327 |
######## END PARAMETERS - CALL FUNCTIONS ##########
|
| 328 |
#[Rates,LumiPageInfo, L1_trig_list,nps]= GetDBRates(run_list, trig_name, trig_list, num_ls, max_dt, physics_active_psi, JSON, debug_print, force_new, SubSystemOff,NoVersion,all_triggers, DoL1,UsePSCol,L1SeedChangeFit)
|
| 329 |
[Rates, LumiPageInfo, L1_trig_list, nps]= GetDBRates(run_list, trig_name, trig_list, num_ls, max_dt, physics_active_psi, JSON, debug_print, force_new, SubSystemOff, NoVersion, all_triggers, DoL1,UsePSCol,L1SeedChangeFit, save_fits)
|
| 330 |
if DoL1:
|
| 331 |
trig_list=L1_trig_list
|
| 332 |
|
| 333 |
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,NoVersion, linear, do_inst,wp_bool,all_triggers,L1SeedChangeFit,nps)
|
| 334 |
|
| 335 |
except KeyboardInterrupt:
|
| 336 |
print "Wait... come back..."
|
| 337 |
|
| 338 |
|
| 339 |
#def GetDBRates(run_list,trig_name,trig_list, num_ls, max_dt, physics_active_psi,JSON,debug_print, force_new, SubSystemOff,NoVersion,all_triggers, DoL1,UsePSCol,L1SeedChangeFit):
|
| 340 |
def GetDBRates(run_list, trig_name, trig_list, num_ls, max_dt, physics_active_psi, JSON, debug_print, force_new, SubSystemOff, NoVersion, all_triggers, DoL1,UsePSCol, L1SeedChangeFit, save_fits):
|
| 341 |
|
| 342 |
Rates = {}
|
| 343 |
LumiPageInfo={}
|
| 344 |
## Save in RefRuns with name dependent on trig_name, num_ls, JSON, and physics_active_psi
|
| 345 |
if JSON:
|
| 346 |
#print "Using JSON file"
|
| 347 |
if physics_active_psi:
|
| 348 |
RefRunNameTemplate = "RefRuns/%s/Rates_%s_%sLS_JPAP.pkl"
|
| 349 |
else:
|
| 350 |
RefRunNameTemplate = "RefRuns/%s/Rates_%s_%sLS_JSON.pkl"
|
| 351 |
else:
|
| 352 |
print "Using Physics and Active ==1"
|
| 353 |
if physics_active_psi:
|
| 354 |
RefRunNameTemplate = "RefRuns/%s/Rates_%s_%sLS_PAP.pkl"
|
| 355 |
else:
|
| 356 |
RefRunNameTemplate = "RefRuns/%s/Rates_%s_%sLS.pkl"
|
| 357 |
|
| 358 |
RefRunFile = RefRunNameTemplate % (thisyear,trig_name,num_ls)
|
| 359 |
RefRunFileHLT = RefRunNameTemplate % (thisyear,"HLT",num_ls)
|
| 360 |
|
| 361 |
print "RefRun=",RefRunFile
|
| 362 |
print "RefRunFileHLT",RefRunFileHLT
|
| 363 |
if not force_new:
|
| 364 |
try: ##Open an existing RefRun file with the same parameters and trigger name
|
| 365 |
pkl_file = open(RefRunFile, 'rb')
|
| 366 |
Rates = pickle.load(pkl_file)
|
| 367 |
pkl_file.close()
|
| 368 |
os.remove(RefRunFile)
|
| 369 |
print "using",RefRunFile
|
| 370 |
|
| 371 |
except:
|
| 372 |
try: ##Open an existing RefRun file with the same parameters and HLT for trigger name
|
| 373 |
pkl_file = open(RefRunFileHLT)
|
| 374 |
HLTRates = pickle.load(pkl_file)
|
| 375 |
for key in HLTRates:
|
| 376 |
if trig_name in str(key):
|
| 377 |
Rates[key] = HLTRates[key]
|
| 378 |
#print str(RefRunFile)+" does not exist. Creating ..."
|
| 379 |
except:
|
| 380 |
print str(RefRunFile)+" does not exist. Creating ..."
|
| 381 |
|
| 382 |
## try the lumis file
|
| 383 |
RefLumiNameTemplate = "RefRuns/%s/Lumis_%s_%sLS.pkl"
|
| 384 |
RefLumiFile= RefLumiNameTemplate % (thisyear,"HLT",num_ls)
|
| 385 |
if not force_new:
|
| 386 |
try:
|
| 387 |
pkl_lumi_file = open(RefLumiFile, 'rb')
|
| 388 |
LumiPageInfo = pickle.load(pkl_lumi_file)
|
| 389 |
pkl_lumi_file.close()
|
| 390 |
os.remove(RefLumiFile)
|
| 391 |
print "using",RefLumiFile
|
| 392 |
except:
|
| 393 |
print str(RefLumiFile)+" doesn't exist. Make it..."
|
| 394 |
|
| 395 |
trig_list_noV=[]
|
| 396 |
for trigs in trig_list:
|
| 397 |
trig_list_noV.append(StripVersion(trigs))
|
| 398 |
|
| 399 |
if NoVersion:
|
| 400 |
trig_list=trig_list_noV
|
| 401 |
|
| 402 |
for RefRunNum in run_list:
|
| 403 |
if JSON:
|
| 404 |
if not RefRunNum in JSON:
|
| 405 |
continue
|
| 406 |
try:
|
| 407 |
ExistsAlready = False
|
| 408 |
for key in Rates:
|
| 409 |
if RefRunNum in Rates[key]["run"]:
|
| 410 |
ExistsAlready = True
|
| 411 |
break
|
| 412 |
|
| 413 |
LumiExistsLAready=False
|
| 414 |
for v in LumiPageInfo.itervalues():
|
| 415 |
if RefRunNum == v["Run"]:
|
| 416 |
LumiExistsAlready=True
|
| 417 |
break
|
| 418 |
if ExistsAlready and LumiExistsAlready:
|
| 419 |
continue
|
| 420 |
|
| 421 |
except:
|
| 422 |
print "Getting info for run "+str(RefRunNum)
|
| 423 |
|
| 424 |
if RefRunNum < 1:
|
| 425 |
continue
|
| 426 |
ColRunNum,isCol,isGood = GetLatestRunNumber(RefRunNum)
|
| 427 |
if not isGood:
|
| 428 |
print "Run ",RefRunNum, " is not Collisions"
|
| 429 |
|
| 430 |
continue
|
| 431 |
|
| 432 |
if not isCol:
|
| 433 |
print "Run ",RefRunNum, " is not Collisions"
|
| 434 |
|
| 435 |
continue
|
| 436 |
|
| 437 |
print "calculating rates and green lumis for run ",RefRunNum
|
| 438 |
|
| 439 |
if True: ##Placeholder
|
| 440 |
if True: #May replace with "try" - for now it's good to know when problems happen
|
| 441 |
RefParser = DatabaseParser()
|
| 442 |
RefParser.RunNumber = RefRunNum
|
| 443 |
RefParser.ParseRunSetup()
|
| 444 |
RefLumiRangePhysicsActive = RefParser.GetLSRange(1,9999) ##Gets array of all LS with physics and active on
|
| 445 |
RefLumiArray = RefParser.GetLumiInfo() ##Gets array of all existing LS and their lumi info
|
| 446 |
RefLumiRange = []
|
| 447 |
RefMoreLumiArray = RefParser.GetMoreLumiInfo()#dict with keys as bits from lumisections WBM page and values are dicts with key=LS:value=bit
|
| 448 |
L1HLTseeds=RefParser.GetL1HLTseeds()
|
| 449 |
HLTL1PS=RefParser.GetL1PSbyseed()
|
| 450 |
###Add all triggers to list if all trigger
|
| 451 |
try:
|
| 452 |
TriggerRatesCheck = RefParser.GetHLTRates([1])##just grab from 1st LS
|
| 453 |
except:
|
| 454 |
print "ERROR: unable to get HLT triggers for this run"
|
| 455 |
exit(2)
|
| 456 |
for HLTkey in TriggerRatesCheck:
|
| 457 |
if NoVersion:
|
| 458 |
name = StripVersion(HLTkey)
|
| 459 |
else:
|
| 460 |
name=HLTkey
|
| 461 |
if not name in trig_list:
|
| 462 |
if all_triggers:
|
| 463 |
trig_list.append(name)
|
| 464 |
|
| 465 |
###add L1 triggers to list if Do L1
|
| 466 |
if DoL1:
|
| 467 |
for HLTkey in trig_list:
|
| 468 |
#print name
|
| 469 |
# if "L1" in HLTkey:
|
| 470 |
# continue
|
| 471 |
if not HLTkey.startswith('HLT'):
|
| 472 |
continue
|
| 473 |
else:
|
| 474 |
try:
|
| 475 |
for L1seed in L1HLTseeds[HLTkey]:
|
| 476 |
if L1seed not in trig_list:
|
| 477 |
trig_list.append(L1seed)
|
| 478 |
except:
|
| 479 |
print "Failed on trigger "+str(HLTkey)
|
| 480 |
pass
|
| 481 |
for iterator in RefLumiArray[0]: ##Makes array of LS with proper PAP and JSON properties
|
| 482 |
##cheap way of getting PSCol None-->0
|
| 483 |
if RefLumiArray[0][iterator] not in range(1,9):
|
| 484 |
RefLumiArray[0][iterator]=0
|
| 485 |
|
| 486 |
if not UsePSCol==-1:
|
| 487 |
if not RefLumiArray[0][iterator]==UsePSCol:
|
| 488 |
print "skipping LS",iterator
|
| 489 |
continue
|
| 490 |
|
| 491 |
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):
|
| 492 |
if not JSON or RefRunNum in JSON:
|
| 493 |
if not JSON or iterator in JSON[RefRunNum]:
|
| 494 |
RefLumiRange.append(iterator)
|
| 495 |
|
| 496 |
try:
|
| 497 |
nls = RefLumiRange[0]
|
| 498 |
LSRange = {}
|
| 499 |
except:
|
| 500 |
print "Run "+str(RefRunNum)+" has no good LS"
|
| 501 |
continue
|
| 502 |
if num_ls > len(RefLumiRange):
|
| 503 |
print "Run "+str(RefRunNum)+" is too short: from "+str(nls)+" to "+str(RefLumiRange[-1])+", while num_ls = "+str(num_ls)
|
| 504 |
continue
|
| 505 |
while nls < RefLumiRange[-1]-num_ls:
|
| 506 |
LSRange[nls] = []
|
| 507 |
counter = 0
|
| 508 |
for iterator in RefLumiRange:
|
| 509 |
if iterator >= nls and counter < num_ls:
|
| 510 |
LSRange[nls].append(iterator)
|
| 511 |
counter += 1
|
| 512 |
nls = LSRange[nls][-1]+1
|
| 513 |
[HLTL1_seedchanges,nps]=checkL1seedChangeALLPScols(trig_list,HLTL1PS) #for L1prescale changes
|
| 514 |
|
| 515 |
|
| 516 |
#print HLTL1_seedchanges
|
| 517 |
#print "nps=",nps
|
| 518 |
#print "Run "+str(RefRunNum)+" contains LS from "+str(min(LSRange))+" to "+str(max(LSRange))
|
| 519 |
for nls in sorted(LSRange.iterkeys()):
|
| 520 |
TriggerRates = RefParser.GetHLTRates(LSRange[nls])
|
| 521 |
#L1Rate=RefParser.GetDeadTimeBeamActive(LSRange[nls])
|
| 522 |
|
| 523 |
## Clumsy way to append Stream A. Should choose correct method for calculating stream a based on ps column used in data taking.
|
| 524 |
|
| 525 |
if ('HLT_Stream_A' in trig_list) or all_triggers:
|
| 526 |
config = RateMonConfig(os.path.abspath(os.path.dirname(sys.argv[0])))
|
| 527 |
config.ReadCFG()
|
| 528 |
stream_mon = StreamMonitor()
|
| 529 |
core_a_rates = stream_mon.getStreamACoreRatesByLS(RefParser,LSRange[nls],config).values()
|
| 530 |
avg_core_a_rate = sum(core_a_rates)/len(LSRange[nls])
|
| 531 |
TriggerRates['HLT_Stream_A'] = [1,1,avg_core_a_rate,avg_core_a_rate]
|
| 532 |
HLTL1_seedchanges["HLT_Stream_A"] = [[ps_col] for ps_col in range(0,nps)]
|
| 533 |
# dummylist=[]
|
| 534 |
# for pscol in range(0,nps):
|
| 535 |
# doubledummylist=[]
|
| 536 |
# doubledummylist.append(pscol)
|
| 537 |
# dummylist.append(doubledummylist)
|
| 538 |
# HLTL1_seedchanges["HLT_Stream_A"]=dummylist
|
| 539 |
|
| 540 |
if DoL1:
|
| 541 |
L1RatesALL=RefParser.GetL1RatesALL(LSRange[nls])
|
| 542 |
|
| 543 |
for L1seed in L1RatesALL.iterkeys():
|
| 544 |
TriggerRates[L1seed]=L1RatesALL[L1seed]
|
| 545 |
|
| 546 |
[inst, live, delivered, dead, pscols] = RefParser.GetAvLumiInfo(LSRange[nls])
|
| 547 |
deadtimebeamactive=RefParser.GetDeadTimeBeamActive(LSRange[nls])
|
| 548 |
|
| 549 |
physics = 1
|
| 550 |
active = 1
|
| 551 |
psi = 99
|
| 552 |
if save_fits and (max(pscols) != min(pscols)):#kick out points which average over two ps columns if doing running in fit making mode
|
| 553 |
continue
|
| 554 |
for iterator in LSRange[nls]: ##Gets lowest value of physics, active, and psi in the set of lumisections
|
| 555 |
if RefLumiArray[5][iterator] == 0:
|
| 556 |
physics = 0
|
| 557 |
if RefLumiArray[6][iterator] == 0:
|
| 558 |
active = 0
|
| 559 |
if RefLumiArray[0][iterator] < psi:
|
| 560 |
psi = RefLumiArray[0][iterator]
|
| 561 |
|
| 562 |
if inst < 0 or live < 0 or delivered < 0:
|
| 563 |
print "Run "+str(RefRunNum)+" LS "+str(nls)+" inst lumi = "+str(inst)+" live lumi = "+str(live)+", delivered = "+str(delivered)+", physics = "+str(physics)+", active = "+str(active)
|
| 564 |
|
| 565 |
LumiPageInfo[nls] = LumiRangeGreens(RefMoreLumiArray,LSRange,nls,RefRunNum,deadtimebeamactive)
|
| 566 |
|
| 567 |
for key in TriggerRates:
|
| 568 |
|
| 569 |
if NoVersion:
|
| 570 |
name = StripVersion(key)
|
| 571 |
else:
|
| 572 |
name=key
|
| 573 |
|
| 574 |
if not name in trig_list:
|
| 575 |
if all_triggers and name.startswith('HLT_Stream_A'):
|
| 576 |
trig_list.append(name) ##Only triggers in trig_list have HLTL1_seedchanges filled
|
| 577 |
else:
|
| 578 |
continue
|
| 579 |
|
| 580 |
if not Rates.has_key(name):
|
| 581 |
Rates[name] = {}
|
| 582 |
Rates[name]["run"] = []
|
| 583 |
Rates[name]["ls"] = []
|
| 584 |
Rates[name]["ps"] = []
|
| 585 |
Rates[name]["inst_lumi"] = []
|
| 586 |
Rates[name]["live_lumi"] = []
|
| 587 |
Rates[name]["delivered_lumi"] = []
|
| 588 |
Rates[name]["deadtime"] = []
|
| 589 |
Rates[name]["rawrate"] = []
|
| 590 |
Rates[name]["rate"] = []
|
| 591 |
Rates[name]["rawxsec"] = []
|
| 592 |
Rates[name]["xsec"] = []
|
| 593 |
Rates[name]["physics"] = []
|
| 594 |
Rates[name]["active"] = []
|
| 595 |
Rates[name]["psi"] = []
|
| 596 |
Rates[name]["L1seedchange"]=[]
|
| 597 |
[avps, ps, rate, psrate] = TriggerRates[key]
|
| 598 |
Rates[name]["run"].append(RefRunNum)
|
| 599 |
Rates[name]["ls"].append(nls)
|
| 600 |
Rates[name]["ps"].append(ps)
|
| 601 |
Rates[name]["inst_lumi"].append(inst)
|
| 602 |
Rates[name]["live_lumi"].append(live)
|
| 603 |
Rates[name]["delivered_lumi"].append(delivered)
|
| 604 |
Rates[name]["deadtime"].append(deadtimebeamactive)
|
| 605 |
Rates[name]["rawrate"].append(rate)
|
| 606 |
Rates[name]["L1seedchange"].append(HLTL1_seedchanges[name])
|
| 607 |
if live == 0:
|
| 608 |
Rates[name]["rate"].append(0.0)
|
| 609 |
Rates[name]["rawxsec"].append(0.0)
|
| 610 |
Rates[name]["xsec"].append(0.0)
|
| 611 |
else:
|
| 612 |
try:
|
| 613 |
Rates[name]["rate"].append(psrate/(1.0-deadtimebeamactive))
|
| 614 |
except:
|
| 615 |
Rates[name]["rate"].append(0.0)
|
| 616 |
Rates[name]["rawxsec"].append(rate/live)
|
| 617 |
Rates[name]["xsec"].append(psrate/live)
|
| 618 |
Rates[name]["physics"].append(physics)
|
| 619 |
Rates[name]["active"].append(active)
|
| 620 |
Rates[name]["psi"].append(psi)
|
| 621 |
|
| 622 |
RateOutput = open(RefRunFile, 'wb') ##Save new Rates[] to RefRuns
|
| 623 |
pickle.dump(Rates, RateOutput, 2)
|
| 624 |
RateOutput.close()
|
| 625 |
LumiOutput = open(RefLumiFile,'wb')
|
| 626 |
pickle.dump(LumiPageInfo,LumiOutput, 2)
|
| 627 |
LumiOutput.close()
|
| 628 |
|
| 629 |
return [Rates,LumiPageInfo,trig_list,nps]
|
| 630 |
|
| 631 |
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,NoVersion, linear,do_inst,wp_bool,all_triggers,L1SeedChangeFit,nps):
|
| 632 |
|
| 633 |
[min_run, max_run, priot, InputFit, OutputFit, OutputFitPS, failed_paths, first_trigger, varX, varY, do_fit, save_root, save_png, fit_file, RootNameTemplate, RootFile, InputFitPS]=InitMakePlots(run_list, trig_name, num_ls, plot_properties, nps, L1SeedChangeFit)
|
| 634 |
##modify for No Version and check the trigger list
|
| 635 |
trig_list=InitTrigList(trig_list, save_fits, NoVersion, InputFit)
|
| 636 |
|
| 637 |
for print_trigger in sorted(Rates):
|
| 638 |
[trig_list, passchecktriglist, meanrawrate] = CheckTrigList(trig_list, print_trigger, all_triggers, masked_triggers, min_rate, Rates, run_list, trig_name, failed_paths)
|
| 639 |
if not passchecktriglist: #failed_paths is modified by CheckTrigList to include output messages explaining why a trigger failed
|
| 640 |
continue
|
| 641 |
|
| 642 |
[meanrate, meanxsec, meanlumi, sloperate, slopexsec, nlow, nhigh, lowrate, lowxsec, lowlumi, highrate, highxsec, highlumi]=GetMeanRates(Rates, print_trigger, max_dt)
|
| 643 |
chioffset=1.0 ##chioffset now a fraction; must be 10% better to use expo rather than quad, quad rather than line
|
| 644 |
width = max([len(trigger_name) for trigger_name in trig_list])
|
| 645 |
for psi in range(0,nps):
|
| 646 |
OutputFitPS[psi][print_trigger]=[]##define empty list for each trigger
|
| 647 |
|
| 648 |
####START OF L1 SEED LOOP####
|
| 649 |
#print "LIST L1 seed changes",Rates[print_trigger]["L1seedchange"][0]
|
| 650 |
|
| 651 |
if L1SeedChangeFit and do_fit:
|
| 652 |
dummyPSColslist=Rates[print_trigger]["L1seedchange"][0]
|
| 653 |
#print print_trigger, dummyPSColslist
|
| 654 |
if len(dummyPSColslist)!=1:
|
| 655 |
dummyPSColslist.append(range(0,nps))
|
| 656 |
else:
|
| 657 |
dummyPSColslist=[]
|
| 658 |
dummyPSColslist.append(range(0,nps))
|
| 659 |
|
| 660 |
|
| 661 |
if not do_fit:
|
| 662 |
[fitparams, passedGetFit, failed_paths, fitparamsPS]=GetFit(do_fit, InputFit, failed_paths, print_trigger, num_ls,L1SeedChangeFit, InputFitPS,nps)
|
| 663 |
if not passedGetFit:
|
| 664 |
print str(print_trigger)+" did not passedGetFit"
|
| 665 |
continue
|
| 666 |
else:
|
| 667 |
fitparams=["unset",0,0,0,0,0,0]
|
| 668 |
fitparamsPS=["unset",{},{},{},{},{},{}]
|
| 669 |
|
| 670 |
for PSColslist in dummyPSColslist:
|
| 671 |
#print print_trigger, PSColslist
|
| 672 |
passPSinCol=0
|
| 673 |
for iterator in range (len(Rates[print_trigger]["run"])):
|
| 674 |
if Rates[print_trigger]["psi"][iterator] in PSColslist:
|
| 675 |
passPSinCol=1
|
| 676 |
#print PSColslist, Rates[print_trigger]["run"][iterator], Rates[print_trigger]["psi"][iterator]
|
| 677 |
if not passPSinCol:
|
| 678 |
##for when there are no LS in some PS col (pretty common!)
|
| 679 |
#print print_trigger, "No data for",PSColslist
|
| 680 |
continue
|
| 681 |
|
| 682 |
|
| 683 |
AllPlotArrays=DoAllPlotArrays(Rates, print_trigger, run_list, data_clean, meanxsec, num_ls, LumiPageInfo, SubSystemOff, max_dt, print_info, trig_list, do_fit, do_inst, debug_print, fitparams, fitparamsPS, L1SeedChangeFit, PSColslist, first_trigger)
|
| 684 |
[VX, VXE, x_label, VY, VYE, y_label, VF, VFE] = GetVXVY(plot_properties, fit_file, AllPlotArrays, L1SeedChangeFit)
|
| 685 |
|
| 686 |
|
| 687 |
####defines gr1 and failure if no graph in OutputFit ####
|
| 688 |
defgrapass = False
|
| 689 |
if len(VX) > 0:
|
| 690 |
[OutputFit,gr1, gr3, failed_paths, defgrapass]=DefineGraphs(print_trigger,OutputFit,do_fit,varX,varY,x_label,y_label,VX,VY,VXE,VYE,VF,VFE,fit_file, failed_paths,PSColslist)
|
| 691 |
if not defgrapass:
|
| 692 |
continue
|
| 693 |
if do_fit:
|
| 694 |
[f1a,f1b,f1c,f1d,first_trigger]= Fitter(gr1,VX,VY,sloperate,nlow,Rates,print_trigger, first_trigger, varX, varY,lowrate)
|
| 695 |
|
| 696 |
|
| 697 |
if print_table or save_fits:
|
| 698 |
###aditional info from f1 params
|
| 699 |
[f1a_Chi2, f1b_Chi2, f1c_Chi2,f1d_Chi2, f1a_BadMinimum, f1b_BadMinimum, f1c_BadMinimum, meanps, av_rte, passmorefitinfo]=more_fit_info(f1a,f1b,f1c,f1d,VX,VY,print_trigger,Rates)
|
| 700 |
if not passmorefitinfo:
|
| 701 |
OutputFit[print_trigger] = ["fit failed","Zero NDF"]
|
| 702 |
###output fit params
|
| 703 |
else:
|
| 704 |
[OutputFit,first_trigger, failed_paths]=output_fit_info(do_fit,f1a,f1b,f1c,f1d,varX,varY,VX,VY,linear,print_trigger,first_trigger,Rates,width,chioffset,wp_bool,num_ls,meanrawrate,OutputFit, failed_paths, PSColslist, dummyPSColslist)
|
| 705 |
if do_fit:
|
| 706 |
for PSI in PSColslist:
|
| 707 |
if not OutputFitPS[PSI][print_trigger]:
|
| 708 |
OutputFitPS[PSI][print_trigger]=OutputFit[print_trigger]
|
| 709 |
|
| 710 |
PSlist=deque(PSColslist)
|
| 711 |
PSmin=PSlist.popleft()
|
| 712 |
if not PSlist:
|
| 713 |
PSmax=PSmin
|
| 714 |
else:
|
| 715 |
PSmax=PSlist.pop()
|
| 716 |
|
| 717 |
first_trigger=False
|
| 718 |
if save_root or save_png:
|
| 719 |
c1 = TCanvas(str(varX),str(varY))
|
| 720 |
c1.SetName(str(print_trigger)+"_ps"+str(PSmin)+"_"+str(PSmax)+"_"+str(varY)+"_vs_"+str(varX))
|
| 721 |
gr1.Draw("APZ")
|
| 722 |
if not do_fit:
|
| 723 |
gr3.Draw("P3")
|
| 724 |
c1.Update()
|
| 725 |
else:
|
| 726 |
c1=DrawFittedCurve(f1a, f1b,f1c, f1d, chioffset,do_fit,c1,VX ,VY,print_trigger,Rates)
|
| 727 |
|
| 728 |
if save_root:
|
| 729 |
myfile = TFile( RootFile, 'UPDATE' )
|
| 730 |
c1.Write()
|
| 731 |
|
| 732 |
|
| 733 |
myfile.Close()
|
| 734 |
if save_png:
|
| 735 |
c1.SaveAs(str(print_trigger)+"_"+str(varY)+"_vs_"+str(varX)+".png")
|
| 736 |
|
| 737 |
|
| 738 |
|
| 739 |
EndMkrootfile(failed_paths, save_fits, save_root, fit_file, RootFile, OutputFit, OutputFitPS, L1SeedChangeFit)
|
| 740 |
|
| 741 |
|
| 742 |
|
| 743 |
|
| 744 |
|
| 745 |
############# SUPPORTING FUNCTIONS ################
|
| 746 |
|
| 747 |
def InitMakePlots(run_list, trig_name, num_ls, plot_properties, nps, L1SeedChangeFit):
|
| 748 |
min_run = min(run_list)
|
| 749 |
max_run = max(run_list)
|
| 750 |
|
| 751 |
priot.has_been_called=False
|
| 752 |
|
| 753 |
InputFit = {}
|
| 754 |
InputFitPS = {}
|
| 755 |
OutputFit = {}
|
| 756 |
failed_paths = []
|
| 757 |
first_trigger=True
|
| 758 |
OutputFitPS={}
|
| 759 |
for ii in range(0,nps):
|
| 760 |
OutputFitPS[ii]={}
|
| 761 |
|
| 762 |
[[varX, varY, do_fit, save_root, save_png, fit_file]] = plot_properties
|
| 763 |
|
| 764 |
RootNameTemplate = "%s_%sLS_%s_vs_%s_Run%s-%s.root"
|
| 765 |
RootFile = RootNameTemplate % (trig_name, num_ls, varX, varY, min_run, max_run)
|
| 766 |
|
| 767 |
if not do_fit:
|
| 768 |
try:
|
| 769 |
pkl_file = open(fit_file, 'rb')
|
| 770 |
InputFit = pickle.load(pkl_file)
|
| 771 |
print "opening fit_file"
|
| 772 |
pkl_file.close()
|
| 773 |
except:
|
| 774 |
print "ERROR: could not open fit file: %s" % (fit_file,)
|
| 775 |
exit(2)
|
| 776 |
if L1SeedChangeFit:
|
| 777 |
try:
|
| 778 |
PSfitfile=fit_file.replace("HLT_NoV","HLT_NoV_ByPS")
|
| 779 |
print "opening",PSfitfile
|
| 780 |
pklfilePS = open(PSfitfile, 'rb')
|
| 781 |
InputFitPS = pickle.load(pklfilePS)
|
| 782 |
except:
|
| 783 |
print "ERROR: could not open fit file: %s" % (PSfitfile,)
|
| 784 |
exit(2)
|
| 785 |
|
| 786 |
if save_root:
|
| 787 |
try:
|
| 788 |
os.remove(RootFile)
|
| 789 |
except:
|
| 790 |
pass
|
| 791 |
|
| 792 |
|
| 793 |
return [min_run, max_run, priot, InputFit, OutputFit, OutputFitPS, failed_paths, first_trigger, varX, varY, do_fit, save_root, save_png, fit_file, RootNameTemplate, RootFile, InputFitPS]
|
| 794 |
|
| 795 |
def InitTrigList(trig_list, save_fits, NoVersion, InputFit):
|
| 796 |
trig_list_noV=[]
|
| 797 |
for trigs in trig_list:
|
| 798 |
trig_list_noV.append(StripVersion(trigs))
|
| 799 |
if NoVersion:
|
| 800 |
trig_list=trig_list_noV
|
| 801 |
|
| 802 |
## check that all the triggers we ask to plot are in the input fit
|
| 803 |
if not save_fits:
|
| 804 |
goodtrig_list = []
|
| 805 |
FitInputNoV={}
|
| 806 |
for trig in trig_list:
|
| 807 |
|
| 808 |
if NoVersion:
|
| 809 |
for trigger in InputFit.iterkeys():
|
| 810 |
FitInputNoV[StripVersion(trigger)]=InputFit[trigger]
|
| 811 |
InputFit=FitInputNoV
|
| 812 |
|
| 813 |
else:
|
| 814 |
if not InputFit.has_key(trig):
|
| 815 |
print "WARNING: No Fit Prediction for Trigger %s, SKIPPING" % (trig,)
|
| 816 |
else:
|
| 817 |
goodtrig_list.append(trig)
|
| 818 |
trig_list = goodtrig_list
|
| 819 |
return trig_list
|
| 820 |
|
| 821 |
##Limits Rates[] to runs in run_list
|
| 822 |
def CheckTrigList(trig_list, print_trigger, all_triggers, masked_triggers, min_rate, Rates, run_list, trig_name, failed_paths):
|
| 823 |
|
| 824 |
NewTrigger = {}
|
| 825 |
passed = 1 ##to replace continue
|
| 826 |
mean_raw_rate = 0
|
| 827 |
if not print_trigger in trig_list:
|
| 828 |
if all_triggers:
|
| 829 |
trig_list.append(print_trigger)
|
| 830 |
else:
|
| 831 |
failed_paths.append([print_trigger,"The monitorlist did not include these paths"])
|
| 832 |
passed = 0
|
| 833 |
return [trig_list, passed, mean_raw_rate]
|
| 834 |
|
| 835 |
for key in Rates[print_trigger]:
|
| 836 |
NewTrigger[key] = []
|
| 837 |
|
| 838 |
for iterator in range(len(Rates[print_trigger]["run"])):
|
| 839 |
if Rates[print_trigger]["run"][iterator] in run_list:
|
| 840 |
for key in Rates[print_trigger]:
|
| 841 |
NewTrigger[key].append(Rates[print_trigger][key][iterator])
|
| 842 |
|
| 843 |
Rates[print_trigger] = NewTrigger
|
| 844 |
mean_raw_rate = sum(Rates[print_trigger]["rawrate"])/len(Rates[print_trigger]["rawrate"])
|
| 845 |
if mean_raw_rate < min_rate:
|
| 846 |
failed_paths.append([print_trigger,"The rate of these paths did not exceed the minimum"])
|
| 847 |
passed = 0
|
| 848 |
|
| 849 |
masked_trig = False
|
| 850 |
for mask in masked_triggers:
|
| 851 |
if str(mask) in print_trigger:
|
| 852 |
masked_trig = True
|
| 853 |
if masked_trig:
|
| 854 |
failed_paths.append([print_trigger,"These paths were masked"])
|
| 855 |
passed = 0
|
| 856 |
|
| 857 |
return [trig_list, passed, mean_raw_rate]
|
| 858 |
|
| 859 |
|
| 860 |
|
| 861 |
def GetMeanRates(Rates, print_trigger, max_dt):
|
| 862 |
lowlumi = 0
|
| 863 |
meanlumi_init = median(Rates[print_trigger]["live_lumi"])
|
| 864 |
meanlumi = 0
|
| 865 |
highlumi = 0
|
| 866 |
lowrate = 0
|
| 867 |
meanrate = 0
|
| 868 |
highrate = 0
|
| 869 |
lowxsec = 0
|
| 870 |
meanxsec = 0
|
| 871 |
highxsec = 0
|
| 872 |
nlow = 0
|
| 873 |
nhigh = 0
|
| 874 |
|
| 875 |
for iterator in range(len(Rates[print_trigger]["rate"])):
|
| 876 |
if Rates[print_trigger]["live_lumi"][iterator] <= meanlumi_init:
|
| 877 |
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):
|
| 878 |
meanrate+=Rates[print_trigger]["rate"][iterator]
|
| 879 |
lowrate+=Rates[print_trigger]["rate"][iterator]
|
| 880 |
meanxsec+=Rates[print_trigger]["xsec"][iterator]
|
| 881 |
lowxsec+=Rates[print_trigger]["xsec"][iterator]
|
| 882 |
meanlumi+=Rates[print_trigger]["live_lumi"][iterator]
|
| 883 |
lowlumi+=Rates[print_trigger]["live_lumi"][iterator]
|
| 884 |
nlow+=1
|
| 885 |
if Rates[print_trigger]["live_lumi"][iterator] > meanlumi_init:
|
| 886 |
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):
|
| 887 |
meanrate+=Rates[print_trigger]["rate"][iterator]
|
| 888 |
highrate+=Rates[print_trigger]["rate"][iterator]
|
| 889 |
meanxsec+=Rates[print_trigger]["xsec"][iterator]
|
| 890 |
highxsec+=Rates[print_trigger]["xsec"][iterator]
|
| 891 |
meanlumi+=Rates[print_trigger]["live_lumi"][iterator]
|
| 892 |
highlumi+=Rates[print_trigger]["live_lumi"][iterator]
|
| 893 |
nhigh+=1
|
| 894 |
try:
|
| 895 |
meanrate = meanrate/(nlow+nhigh)
|
| 896 |
meanxsec = meanxsec/(nlow+nhigh)
|
| 897 |
meanlumi = meanlumi/(nlow+nhigh)
|
| 898 |
if (nlow==0):
|
| 899 |
sloperate = (highrate/nhigh) / (highlumi/nhigh)
|
| 900 |
slopexsec = (highxsec/nhigh) / (highlumi/nhigh)
|
| 901 |
elif (nhigh==0):
|
| 902 |
sloperate = (lowrate/nlow) / (lowlumi/nlow)
|
| 903 |
slopexsec = (lowxsec/nlow) / (lowlumi/nlow)
|
| 904 |
else:
|
| 905 |
sloperate = ( (highrate/nhigh) - (lowrate/nlow) ) / ( (highlumi/nhigh) - (lowlumi/nlow) )
|
| 906 |
slopexsec = ( (highxsec/nhigh) - (lowxsec/nlow) ) / ( (highlumi/nhigh) - (lowlumi/nlow) )
|
| 907 |
except:
|
| 908 |
# print str(print_trigger)+" has no good datapoints - setting initial xsec slope estimate to 0"
|
| 909 |
meanrate = median(Rates[print_trigger]["rate"])
|
| 910 |
meanxsec = median(Rates[print_trigger]["xsec"])
|
| 911 |
meanlumi = median(Rates[print_trigger]["live_lumi"])
|
| 912 |
sloperate = meanxsec
|
| 913 |
slopexsec = 0
|
| 914 |
return [meanrate, meanxsec, meanlumi, sloperate, slopexsec, nlow, nhigh, lowrate, lowxsec, lowlumi, highrate, highxsec, highlumi]
|
| 915 |
|
| 916 |
|
| 917 |
################
|
| 918 |
def GetFit(do_fit, InputFit, failed_paths, print_trigger, num_ls, L1SeedChangeFit, InputFitPS,nps):
|
| 919 |
|
| 920 |
passed=1
|
| 921 |
FitTypePS={}
|
| 922 |
X0PS={}
|
| 923 |
X1PS={}
|
| 924 |
X2PS={}
|
| 925 |
X3PS={}
|
| 926 |
sigmaPS={}
|
| 927 |
X0errPS={}
|
| 928 |
|
| 929 |
try:
|
| 930 |
FitType = InputFit[print_trigger][0]
|
| 931 |
except:
|
| 932 |
failed_paths.append([print_trigger,"These paths did not exist in the monitorlist used to create the fit"])
|
| 933 |
FitType = "parse failed"
|
| 934 |
passed=0
|
| 935 |
fitparams=[FitType, 0, 0, 0, 0, 0, 0]
|
| 936 |
if FitType == "fit failed":
|
| 937 |
failure_comment = InputFit[print_trigger][1]
|
| 938 |
failed_paths.append([print_trigger, failure_comment])
|
| 939 |
passed=0
|
| 940 |
fitparams=[FitType, 0, 0, 0, 0, 0, 0]
|
| 941 |
elif FitType == "parse failed":
|
| 942 |
failure_comment = "These paths did not exist in the monitorlist used to create the fit"
|
| 943 |
else:
|
| 944 |
X0 = InputFit[print_trigger][1]
|
| 945 |
X1 = InputFit[print_trigger][2]
|
| 946 |
X2 = InputFit[print_trigger][3]
|
| 947 |
X3 = InputFit[print_trigger][4]
|
| 948 |
sigma = InputFit[print_trigger][5]/math.sqrt(num_ls)*3#Display 3 sigma band to show outliers more clearly
|
| 949 |
X0err= InputFit[print_trigger][7]
|
| 950 |
fitparams=[FitType, X0, X1, X2, X3, sigma, X0err]
|
| 951 |
|
| 952 |
|
| 953 |
if L1SeedChangeFit:
|
| 954 |
|
| 955 |
for psi in range(0,nps):
|
| 956 |
#print psi, print_trigger, InputFitPS[psi][print_trigger]
|
| 957 |
try:
|
| 958 |
FitTypePS[psi] = InputFitPS[psi][print_trigger][0]
|
| 959 |
except:
|
| 960 |
failed_paths.append([print_trigger+'_PS_'+str(psi),"These paths did not exist in the monitorlist used to create the fit"])
|
| 961 |
FitTypePS[psi] = "parse failed"
|
| 962 |
passed=0
|
| 963 |
fitparamsPS=[FitTypePS, X0PS, X1PS, X2PS, X3PS, sigmaPS, X0errPS]
|
| 964 |
if FitTypePS[psi] == "fit failed":
|
| 965 |
failure_comment = InputFitPS[psi][print_trigger][1]
|
| 966 |
failed_paths.append([print_trigger+str(psi), failure_comment])
|
| 967 |
passed=0
|
| 968 |
fitparamsPS=[FitTypePS, X0PS, X1PS, X2PS, X3PS, sigmaPS, X0errPS]
|
| 969 |
else:
|
| 970 |
try:
|
| 971 |
X0PS[psi] = InputFitPS[psi][print_trigger][1]
|
| 972 |
X1PS[psi] = InputFitPS[psi][print_trigger][2]
|
| 973 |
X2PS[psi] = InputFitPS[psi][print_trigger][3]
|
| 974 |
X3PS[psi] = InputFitPS[psi][print_trigger][4]
|
| 975 |
sigmaPS[psi] = InputFitPS[psi][print_trigger][5]/math.sqrt(num_ls)*3#Display 3 sigma band to show outliers more clearly
|
| 976 |
X0errPS[psi]= InputFitPS[psi][print_trigger][7]
|
| 977 |
fitparamsPS=[FitTypePS, X0PS, X1PS, X2PS, X3PS, sigmaPS, X0errPS]
|
| 978 |
except:
|
| 979 |
#print "ERROR: unable to get fits by PS for",print_trigger," in col",psi, "skipping."
|
| 980 |
pass
|
| 981 |
|
| 982 |
|
| 983 |
|
| 984 |
return [fitparams, passed, failed_paths, fitparamsPS]
|
| 985 |
|
| 986 |
## we are 2 lumis off when we start! -gets worse when we skip lumis
|
| 987 |
def DoAllPlotArrays(Rates, print_trigger, run_list, data_clean, meanxsec, num_ls, LumiPageInfo, SubSystemOff, max_dt, print_info, trig_list, do_fit, do_inst, debug_print, fitparams, fitparamsPS, L1SeedChangeFit, PSColslist, first_trigger):
|
| 988 |
|
| 989 |
###init arrays ###
|
| 990 |
[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()
|
| 991 |
|
| 992 |
|
| 993 |
|
| 994 |
it_offset=0
|
| 995 |
###loop over each LS ###
|
| 996 |
for iterator in range(len(Rates[print_trigger]["rate"])):
|
| 997 |
if not Rates[print_trigger]["run"][iterator] in run_list:
|
| 998 |
continue
|
| 999 |
if not Rates[print_trigger]["psi"][iterator] in PSColslist:
|
| 1000 |
continue
|
| 1001 |
|
| 1002 |
#else:
|
| 1003 |
#print iterator, Rates[print_trigger]["psi"][iterator], PSColslist
|
| 1004 |
##Old Spike-killer: used average-xsec method, too clumsy, esp. for nonlinear triggers
|
| 1005 |
#prediction = meanxsec + slopexsec * (Rates[print_trigger]["live_lumi"][iterator] - meanlumi)
|
| 1006 |
#realvalue = Rates[print_trigger]["xsec"][iterator]
|
| 1007 |
|
| 1008 |
##New Spike-killer: gets average of nearest 4 LS
|
| 1009 |
try:
|
| 1010 |
if (iterator > 2 and iterator+2 < len(Rates[print_trigger]["rate"])): ##2 LS before, 2 after
|
| 1011 |
prediction = (Rates[print_trigger]["rate"][iterator-2]+Rates[print_trigger]["rate"][iterator-1]+Rates[print_trigger]["rate"][iterator+1]+Rates[print_trigger]["rate"][iterator+2])/4.0
|
| 1012 |
elif (iterator > 2 and len(Rates[print_trigger]["rate"]) > 4 ): ##4 LS before
|
| 1013 |
prediction = (Rates[print_trigger]["rate"][iterator-4]+Rates[print_trigger]["rate"][iterator-3]+Rates[print_trigger]["rate"][iterator-2]+Rates[print_trigger]["rate"][iterator-1])/4.0
|
| 1014 |
elif (iterator+2 < len(Rates[print_trigger]["rate"]) and len(Rates[print_trigger]["rate"]) > 4 ): ##4 LS after
|
| 1015 |
prediction = (Rates[print_trigger]["rate"][iterator+1]+Rates[print_trigger]["rate"][iterator+2]+Rates[print_trigger]["rate"][iterator+3]+Rates[print_trigger]["rate"][iterator+4])/4.0
|
| 1016 |
else:
|
| 1017 |
prediction = Rates[print_trigger]["rate"][iterator]
|
| 1018 |
realvalue = Rates[print_trigger]["rate"][iterator]
|
| 1019 |
except:
|
| 1020 |
print "Error calculating prediction. Setting rates to defaults."
|
| 1021 |
prediction = Rates[print_trigger]["rate"][iterator]
|
| 1022 |
realvalue = Rates[print_trigger]["rate"][iterator]
|
| 1023 |
|
| 1024 |
if pass_cuts(data_clean, realvalue, prediction, meanxsec, Rates, print_trigger, iterator, num_ls,LumiPageInfo,SubSystemOff,max_dt,print_info, trig_list, first_trigger):
|
| 1025 |
run_t.append(Rates[print_trigger]["run"][iterator])
|
| 1026 |
ls_t.append(Rates[print_trigger]["ls"][iterator])
|
| 1027 |
ps_t.append(Rates[print_trigger]["ps"][iterator])
|
| 1028 |
inst_t.append(Rates[print_trigger]["inst_lumi"][iterator])
|
| 1029 |
live_t.append(Rates[print_trigger]["live_lumi"][iterator])
|
| 1030 |
delivered_t.append(Rates[print_trigger]["delivered_lumi"][iterator])
|
| 1031 |
deadtime_t.append(Rates[print_trigger]["deadtime"][iterator])
|
| 1032 |
rawrate_t.append(Rates[print_trigger]["rawrate"][iterator])
|
| 1033 |
rate_t.append(Rates[print_trigger]["rate"][iterator])
|
| 1034 |
rawxsec_t.append(Rates[print_trigger]["rawxsec"][iterator])
|
| 1035 |
xsec_t.append(Rates[print_trigger]["xsec"][iterator])
|
| 1036 |
psi_t.append(Rates[print_trigger]["psi"][iterator])
|
| 1037 |
|
| 1038 |
e_run_t.append(0.0)
|
| 1039 |
e_ls_t.append(0.0)
|
| 1040 |
e_ps_t.append(0.0)
|
| 1041 |
e_inst_t.append(14.14)
|
| 1042 |
e_live_t.append(14.14)
|
| 1043 |
e_delivered_t.append(14.14)
|
| 1044 |
e_deadtime_t.append(0.01)
|
| 1045 |
e_rawrate_t.append(math.sqrt(Rates[print_trigger]["rawrate"][iterator]/(num_ls*23.3)))
|
| 1046 |
e_rate_t.append(Rates[print_trigger]["ps"][iterator]*math.sqrt(Rates[print_trigger]["rawrate"][iterator]/(num_ls*23.3)))
|
| 1047 |
e_psi_t.append(0.0)
|
| 1048 |
|
| 1049 |
if live_t[-1] == 0:
|
| 1050 |
e_rawxsec_t.append(0)
|
| 1051 |
e_xsec_t.append(0)
|
| 1052 |
else:
|
| 1053 |
try:
|
| 1054 |
e_rawxsec_t.append(math.sqrt(Rates[print_trigger]["rawrate"][iterator]/(num_ls*23.3))/Rates[print_trigger]["live_lumi"][iterator])
|
| 1055 |
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])
|
| 1056 |
except:
|
| 1057 |
e_rawxsec_t.append(0.)
|
| 1058 |
e_xsec_t.append(0.)
|
| 1059 |
|
| 1060 |
if not do_fit:
|
| 1061 |
[FitType, X0, X1, X2, X3, sigma, X0err] = GetCorrectFitParams(fitparams,fitparamsPS,Rates,L1SeedChangeFit,iterator,print_trigger)
|
| 1062 |
if not do_inst:
|
| 1063 |
if FitType == "expo":
|
| 1064 |
rate_prediction = X0 + X1*math.exp(X2+X3*delivered_t[-1])
|
| 1065 |
else:
|
| 1066 |
rate_prediction = X0 + X1*delivered_t[-1] + X2*delivered_t[-1]*delivered_t[-1] + X3*delivered_t[-1]*delivered_t[-1]*delivered_t[-1]
|
| 1067 |
|
| 1068 |
else:
|
| 1069 |
if FitType == "expo":
|
| 1070 |
rate_prediction = X0 + X1*math.exp(X2+X3*inst_t[-1])
|
| 1071 |
else:
|
| 1072 |
rate_prediction = X0 + X1*inst_t[-1] + X2*inst_t[-1]*inst_t[-1] + X3*inst_t[-1]*inst_t[-1]*inst_t[-1]
|
| 1073 |
|
| 1074 |
if rate_prediction != abs(rate_prediction):
|
| 1075 |
rate_prediction = 0
|
| 1076 |
print 'Problem calculating rate prediction. Setting to 0 for '+print_trigger+': lumisection '+ls_t[-1]
|
| 1077 |
|
| 1078 |
if live_t[-1] == 0:
|
| 1079 |
rawrate_fit_t.append(0)
|
| 1080 |
rate_fit_t.append(0)
|
| 1081 |
rawxsec_fit_t.append(0)
|
| 1082 |
xsec_fit_t.append(0)
|
| 1083 |
e_rawrate_fit_t.append(0)
|
| 1084 |
e_rate_fit_t.append(sigma)
|
| 1085 |
e_rawxsec_fit_t.append(0)
|
| 1086 |
e_xsec_fit_t.append(0)
|
| 1087 |
|
| 1088 |
else:
|
| 1089 |
if ps_t[-1]>0.0:
|
| 1090 |
rawrate_fit_t.append(rate_prediction*(1.0-deadtime_t[-1])/(ps_t[-1]))
|
| 1091 |
else:
|
| 1092 |
rawrate_fit_t.append(0.0)
|
| 1093 |
|
| 1094 |
rate_fit_t.append(rate_prediction)
|
| 1095 |
e_rate_fit_t.append(sigma*math.sqrt(rate_prediction))
|
| 1096 |
rawxsec_fit_t.append(rawrate_fit_t[-1]/live_t[-1])
|
| 1097 |
xsec_fit_t.append(rate_prediction*(1.0-deadtime_t[-1])/live_t[-1])
|
| 1098 |
try:
|
| 1099 |
e_rawrate_fit_t.append(sigma*math.sqrt(rate_fit_t[-1])*rawrate_fit_t[-1]/rate_fit_t[-1])
|
| 1100 |
e_rawxsec_fit_t.append(sigma*math.sqrt(rate_fit_t[-1])*rawxsec_fit_t[-1]/rate_fit_t[-1])
|
| 1101 |
e_xsec_fit_t.append(sigma*math.sqrt(rate_fit_t[-1])*xsec_fit_t[-1]/rate_fit_t[-1])
|
| 1102 |
except:
|
| 1103 |
print print_trigger, "has no fitted rate for LS", Rates[print_trigger]["ls"][iterator]
|
| 1104 |
e_rawrate_fit_t.append(sigma)
|
| 1105 |
e_rawxsec_fit_t.append(sigma)
|
| 1106 |
e_xsec_fit_t.append(sigma)
|
| 1107 |
|
| 1108 |
|
| 1109 |
if (print_info and num_ls==1 and (fabs(rawrate_fit_t[-1]-rawrate_t[-1])>2.5*sqrt(sum(Rates[print_trigger]["rawrate"])/len(Rates[print_trigger]["rawrate"])))):
|
| 1110 |
pass
|
| 1111 |
|
| 1112 |
else: ##If the data point does not pass the data_clean filter
|
| 1113 |
if debug_print:
|
| 1114 |
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)
|
| 1115 |
|
| 1116 |
## End "for iterator in range(len(Rates[print_trigger]["rate"])):" loop
|
| 1117 |
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]
|
| 1118 |
|
| 1119 |
def GetCorrectFitParams(fitparams,fitparamsPS,Rates,L1SeedChangeFit,iterator,print_trigger):
|
| 1120 |
if not L1SeedChangeFit:
|
| 1121 |
return fitparams
|
| 1122 |
else:
|
| 1123 |
psi=Rates[print_trigger]["psi"][iterator]
|
| 1124 |
[FitTypePS, X0PS, X1PS, X2PS, X3PS, sigmaPS, X0errPS]=fitparamsPS
|
| 1125 |
return [FitTypePS[psi], X0PS[psi], X1PS[psi], X2PS[psi], X3PS[psi], sigmaPS[psi], X0errPS[psi]]
|
| 1126 |
#[FitType, X0, X1, X2, X3, sigma, X0err]
|
| 1127 |
|
| 1128 |
|
| 1129 |
def CalcSigma(var_x, var_y, func, do_high_lumi):
|
| 1130 |
residuals = []
|
| 1131 |
residuals_high_lumi = []
|
| 1132 |
res_frac = []
|
| 1133 |
res_frac_high_lumi = []
|
| 1134 |
for x, y in zip(var_x,var_y):
|
| 1135 |
y_predicted = func.Eval(x,0,0)
|
| 1136 |
residuals.append(y - y_predicted)
|
| 1137 |
res_frac.append((y - y_predicted)/math.sqrt(abs(y_predicted))) #QUICK FIX, WE NEED TO DECIDE HOW TO HANDLE NEGATIVE
|
| 1138 |
if x > 6000:
|
| 1139 |
residuals_high_lumi.append(y - y_predicted)
|
| 1140 |
res_frac_high_lumi.append((y - y_predicted)/math.sqrt(abs(y_predicted)))
|
| 1141 |
|
| 1142 |
res_squared = [i*i for i in residuals]
|
| 1143 |
res_frac_squared = [i*i for i in res_frac]
|
| 1144 |
res_high_lumi_squared = [i*i for i in residuals_high_lumi]
|
| 1145 |
res_frac_high_lumi_squared = [i*i for i in res_frac_high_lumi]
|
| 1146 |
dev_high_lumi_squared = [i*fabs(i) for i in residuals_high_lumi]
|
| 1147 |
dev_frac_high_lumi_squared = [i*fabs(i) for i in res_frac_high_lumi]
|
| 1148 |
|
| 1149 |
if len(res_squared) > 2:
|
| 1150 |
sigma = math.sqrt(sum(res_squared)/(1.0*len(res_squared)-2.0))
|
| 1151 |
sigma_frac = math.sqrt(sum(res_frac_squared)/(1.0*len(res_frac_squared)-2.0))
|
| 1152 |
else:
|
| 1153 |
sigma = 0
|
| 1154 |
sigma_frac = 0
|
| 1155 |
|
| 1156 |
if len(res_high_lumi_squared) > 10 and do_high_lumi:
|
| 1157 |
high_lumi_sigma_frac = math.sqrt(sum(res_frac_high_lumi_squared)/(1.0*len(res_frac_high_lumi_squared))) ##Statistics limited, don't subtract 2
|
| 1158 |
high_lumi_dev_frac = math.sqrt( fabs( sum(dev_frac_high_lumi_squared)/(1.0*len(dev_frac_high_lumi_squared)) ) ) ##Statistics limited, don't subtract 2
|
| 1159 |
if high_lumi_sigma_frac > 1.25*sigma_frac:
|
| 1160 |
#print "high_lumi_sigma_frac is higher by "+str(100*round((high_lumi_sigma_frac/sigma_frac)-1,2))+"% than sigma_frac ("+str(round(sigma_frac,2))+")"
|
| 1161 |
sigma = sigma*( 0.5 + 0.5*(high_lumi_sigma_frac/sigma_frac) )
|
| 1162 |
sigma_frac = sigma_frac*( 0.5 + 0.5*(high_lumi_sigma_frac/sigma_frac) )
|
| 1163 |
if high_lumi_dev_frac > 4.0*math.sqrt(1.0/(1.0*len(res_frac_high_lumi_squared)-2.0))*sigma_frac:
|
| 1164 |
#print "Total points: "+str(len(res_frac_squared))
|
| 1165 |
#print "High lumi points: "+str(len(res_frac_high_lumi_squared))
|
| 1166 |
#print "high_lumi_dev_frac is "+str(100*round(high_lumi_dev_frac/sigma_frac,2))+"% of sigma_frac ("+str(round(sigma_frac,2))+")"
|
| 1167 |
sigma = sigma*(1.0 + 0.5*(high_lumi_dev_frac/sigma_frac) )
|
| 1168 |
sigma_frac = sigma_frac*(1.0 + 0.5*(high_lumi_dev_frac/sigma_frac) )
|
| 1169 |
|
| 1170 |
return sigma_frac
|
| 1171 |
|
| 1172 |
def GetJSON(json_file):
|
| 1173 |
|
| 1174 |
input_file = open(json_file)
|
| 1175 |
file_content = input_file.read()
|
| 1176 |
inputRange = selectionParser(file_content)
|
| 1177 |
JSON = inputRange.runsandls()
|
| 1178 |
return JSON
|
| 1179 |
##JSON is an array: JSON[run_number] = [1st ls, 2nd ls, 3rd ls ... nth ls]
|
| 1180 |
|
| 1181 |
def MakePlotArrays():
|
| 1182 |
run_t = array.array('f')
|
| 1183 |
ls_t = array.array('f')
|
| 1184 |
ps_t = array.array('f')
|
| 1185 |
inst_t = array.array('f')
|
| 1186 |
live_t = array.array('f')
|
| 1187 |
delivered_t = array.array('f')
|
| 1188 |
deadtime_t = array.array('f')
|
| 1189 |
rawrate_t = array.array('f')
|
| 1190 |
rate_t = array.array('f')
|
| 1191 |
rawxsec_t = array.array('f')
|
| 1192 |
xsec_t = array.array('f')
|
| 1193 |
psi_t = array.array('f')
|
| 1194 |
|
| 1195 |
e_run_t = array.array('f')
|
| 1196 |
e_ls_t = array.array('f')
|
| 1197 |
e_ps_t = array.array('f')
|
| 1198 |
e_inst_t = array.array('f')
|
| 1199 |
e_live_t = array.array('f')
|
| 1200 |
e_delivered_t = array.array('f')
|
| 1201 |
e_deadtime_t = array.array('f')
|
| 1202 |
e_rawrate_t = array.array('f')
|
| 1203 |
e_rate_t = array.array('f')
|
| 1204 |
e_rawxsec_t = array.array('f')
|
| 1205 |
e_xsec_t = array.array('f')
|
| 1206 |
e_psi_t = array.array('f')
|
| 1207 |
|
| 1208 |
rawrate_fit_t = array.array('f')
|
| 1209 |
rate_fit_t = array.array('f')
|
| 1210 |
rawxsec_fit_t = array.array('f')
|
| 1211 |
xsec_fit_t = array.array('f')
|
| 1212 |
e_rawrate_fit_t = array.array('f')
|
| 1213 |
e_rate_fit_t = array.array('f')
|
| 1214 |
e_rawxsec_fit_t = array.array('f')
|
| 1215 |
e_xsec_fit_t = array.array('f')
|
| 1216 |
|
| 1217 |
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]
|
| 1218 |
|
| 1219 |
|
| 1220 |
def GetVXVY(plot_properties, fit_file, AllPlotArrays, L1SeedChangeFit):
|
| 1221 |
|
| 1222 |
VF = "0"
|
| 1223 |
VFE = "0"
|
| 1224 |
|
| 1225 |
[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
|
| 1226 |
for varX, varY, do_fit, save_root, save_png, fit_file in plot_properties:
|
| 1227 |
if varX == "run":
|
| 1228 |
VX = run_t
|
| 1229 |
VXE = run_t_e
|
| 1230 |
x_label = "Run Number"
|
| 1231 |
elif varX == "ls":
|
| 1232 |
VX = ls_t
|
| 1233 |
VXE = e_ls_t
|
| 1234 |
x_label = "Lumisection"
|
| 1235 |
elif varX == "ps":
|
| 1236 |
VX = ps_t
|
| 1237 |
VXE = e_ps_t
|
| 1238 |
x_label = "Prescale"
|
| 1239 |
elif varX == "inst":
|
| 1240 |
VX = inst_t
|
| 1241 |
VXE = e_inst_t
|
| 1242 |
x_label = "Instantaneous Luminosity [10^{30} Hz/cm^{2}]"
|
| 1243 |
|
| 1244 |
elif varX == "live":
|
| 1245 |
VX = live_t
|
| 1246 |
VXE = e_live_t
|
| 1247 |
x_label = "Instantaneous Luminosity [10^{30} Hz/cm^{2}]"
|
| 1248 |
|
| 1249 |
elif varX == "delivered":
|
| 1250 |
VX = delivered_t
|
| 1251 |
VXE = e_delivered_t
|
| 1252 |
x_label = "Instantaneous Luminosity [10^{30} Hz/cm^{2}]"
|
| 1253 |
|
| 1254 |
elif varX == "deadtime":
|
| 1255 |
VX = deadtime_t
|
| 1256 |
VXE = e_deadtime_t
|
| 1257 |
x_label = "Deadtime"
|
| 1258 |
elif varX == "rawrate":
|
| 1259 |
VX = rawrate_t
|
| 1260 |
VXE = e_rawrate_t
|
| 1261 |
x_label = "Raw Rate [Hz]"
|
| 1262 |
elif varX == "rate":
|
| 1263 |
VX = rate_t
|
| 1264 |
VXE = e_rate_t
|
| 1265 |
x_label = "Rate [Hz]"
|
| 1266 |
elif varX == "rawxsec":
|
| 1267 |
VX = rawxsec_t
|
| 1268 |
VXE = e_rawxsec_t
|
| 1269 |
x_label = "Cross Section"
|
| 1270 |
elif varX == "xsec":
|
| 1271 |
VX = xsec_t
|
| 1272 |
VXE = e_xsec_t
|
| 1273 |
x_label = "Cross Section"
|
| 1274 |
elif varX == "psi":
|
| 1275 |
VX = psi_t
|
| 1276 |
VXE = e_psi_t
|
| 1277 |
x_label = "Prescale Index"
|
| 1278 |
else:
|
| 1279 |
print "No valid variable entered for X"
|
| 1280 |
continue
|
| 1281 |
if varY == "run":
|
| 1282 |
VY = run_t
|
| 1283 |
VYE = run_t_e
|
| 1284 |
y_label = "Run Number"
|
| 1285 |
elif varY == "ls":
|
| 1286 |
VY = ls_t
|
| 1287 |
VYE = e_ls_t
|
| 1288 |
y_label = "Lumisection"
|
| 1289 |
elif varY == "ps":
|
| 1290 |
VY = ps_t
|
| 1291 |
VYE = e_ps_t
|
| 1292 |
y_label = "Prescale"
|
| 1293 |
elif varY == "inst":
|
| 1294 |
VY = inst_t
|
| 1295 |
VYE = e_inst_t
|
| 1296 |
y_label = "Instantaneous Luminosity"
|
| 1297 |
elif varY == "live":
|
| 1298 |
VY = live_t
|
| 1299 |
VYE = e_live_t
|
| 1300 |
y_label = "Instantaneous Luminosity"
|
| 1301 |
elif varY == "delivered":
|
| 1302 |
VY = delivered_t
|
| 1303 |
VYE = e_delivered_t
|
| 1304 |
y_label = "Instantaneous Luminosity"
|
| 1305 |
elif varY == "deadtime":
|
| 1306 |
VY = deadtime_t
|
| 1307 |
VYE = e_deadtime_t
|
| 1308 |
y_label = "Deadtime"
|
| 1309 |
elif varY == "rawrate":
|
| 1310 |
VY = rawrate_t
|
| 1311 |
VYE = e_rawrate_t
|
| 1312 |
y_label = "Raw Rate [Hz]"
|
| 1313 |
if fit_file:
|
| 1314 |
VF = rawrate_fit_t
|
| 1315 |
VFE = e_rawrate_fit_t
|
| 1316 |
elif varY == "rate":
|
| 1317 |
VY = rate_t
|
| 1318 |
VYE = e_rate_t
|
| 1319 |
y_label = "Rate [Hz]"
|
| 1320 |
if fit_file:
|
| 1321 |
VF = rate_fit_t
|
| 1322 |
VFE = e_rate_fit_t
|
| 1323 |
elif varY == "rawxsec":
|
| 1324 |
VY = rawxsec_t
|
| 1325 |
VYE = e_rawxsec_t
|
| 1326 |
y_label = "Cross Section"
|
| 1327 |
if fit_file:
|
| 1328 |
VF = rawxsec_fit_t
|
| 1329 |
VFE = e_rawxsec_fit_t
|
| 1330 |
elif varY == "xsec":
|
| 1331 |
VY = xsec_t
|
| 1332 |
VYE = e_xsec_t
|
| 1333 |
y_label = "Cross Section"
|
| 1334 |
if fit_file:
|
| 1335 |
VF = xsec_fit_t
|
| 1336 |
VFE = e_xsec_fit_t
|
| 1337 |
elif varY == "psi":
|
| 1338 |
VY = psi_t
|
| 1339 |
VYE = e_psi_t
|
| 1340 |
y_label = "Prescale Index"
|
| 1341 |
else:
|
| 1342 |
print "No valid variable entered for Y"
|
| 1343 |
continue
|
| 1344 |
|
| 1345 |
return [VX, VXE, x_label, VY, VYE, y_label, VF, VFE]
|
| 1346 |
|
| 1347 |
def pass_cuts(data_clean, realvalue, prediction, meanxsec, Rates, print_trigger, iterator, num_ls,LumiPageInfo,SubSystemOff, max_dt, print_info, trig_list, first_trigger):
|
| 1348 |
it_offset=0
|
| 1349 |
Passed=True
|
| 1350 |
subsystemfailed=[]
|
| 1351 |
|
| 1352 |
if num_ls==1:
|
| 1353 |
##fit is 2 ls ahead of real rate
|
| 1354 |
LS=Rates[print_trigger]["ls"][iterator]
|
| 1355 |
LSRange=LumiPageInfo[LS]["LSRange"]
|
| 1356 |
LS2=LSRange[-1]
|
| 1357 |
lumidict={}
|
| 1358 |
lumidict=LumiPageInfo[LS]
|
| 1359 |
|
| 1360 |
if print_info:
|
| 1361 |
if (iterator==0 and print_trigger==trig_list[0] and first_trigger):
|
| 1362 |
print '%10s%10s%10s%10s%10s%10s%10s%15s%20s%15s' % ("Status", "Run", "LS", "Physics", "Active", "Deadtime", " MaxDeadTime", " Passed all subsystems?", " List of Subsystems", " Spike killing")
|
| 1363 |
|
| 1364 |
## if SubSystemOff["All"]:
|
| 1365 |
## for keys in LumiPageInfo[LS]:
|
| 1366 |
## #print LS, keys, LumiPageInfo[LS][keys]
|
| 1367 |
## if not LumiPageInfo[LS][keys]:
|
| 1368 |
## Passed=False
|
| 1369 |
## subsystemfailed.append(keys)
|
| 1370 |
## break
|
| 1371 |
## else:
|
| 1372 |
if SubSystemOff["Mu"] or SubSystemOff["All"]:
|
| 1373 |
if not (LumiPageInfo[LS]["rpc"] and LumiPageInfo[LS]["dt0"] and LumiPageInfo[LS]["dtp"] and LumiPageInfo[LS]["dtm"] and LumiPageInfo[LS]["cscp"] and LumiPageInfo[LS]["cscm"]):
|
| 1374 |
Passed=False
|
| 1375 |
subsystemfailed.append("Mu")
|
| 1376 |
if SubSystemOff["HCal"] or SubSystemOff["All"]:
|
| 1377 |
if not (LumiPageInfo[LS]["hbhea"] and LumiPageInfo[LS]["hbheb"] and LumiPageInfo[LS]["hbhec"]):
|
| 1378 |
Passed=False
|
| 1379 |
subsystemfailed.append("HCal")
|
| 1380 |
if (SubSystemOff["EndCap"] or SubSystemOff["All"]) and not (LumiPageInfo[LS]["hf"]):
|
| 1381 |
Passed=False
|
| 1382 |
subsystemfailed.append("HCal-EndCap")
|
| 1383 |
if SubSystemOff["ECal"] or SubSystemOff["All"]:
|
| 1384 |
if not (LumiPageInfo[LS]["ebp"] and LumiPageInfo[LS]["ebm"]):
|
| 1385 |
Passed=False
|
| 1386 |
subsystemfailed.append("ECal")
|
| 1387 |
if (SubSystemOff["EndCap"] or SubSystemOff["All"]) and not (LumiPageInfo[LS]["eep"] and LumiPageInfo[LS]["eem"] and LumiPageInfo[LS]["esp"] or LumiPageInfo[LS]["esm"]):
|
| 1388 |
Passed=False
|
| 1389 |
subsystemfailed.append("ECal-EndCap")
|
| 1390 |
if SubSystemOff["Tracker"] or SubSystemOff["All"]:
|
| 1391 |
if not (LumiPageInfo[LS]["tob"] and LumiPageInfo[LS]["tibtid"] and LumiPageInfo[LS]["bpix"] and LumiPageInfo[LS]["fpix"]):
|
| 1392 |
Passed=False
|
| 1393 |
subsystemfailed.append("Tracker")
|
| 1394 |
if (SubSystemOff["EndCap"] or SubSystemOff["All"]) and not (LumiPageInfo[LS]["tecp"] and LumiPageInfo[LS]["tecm"]):
|
| 1395 |
Passed=False
|
| 1396 |
subsystemfailed.append("Tracker-EndCap")
|
| 1397 |
if SubSystemOff["Beam"] or SubSystemOff["All"]:
|
| 1398 |
if not(LumiPageInfo[LS]["b1pres"] and LumiPageInfo[LS]["b2pres"] and LumiPageInfo[LS]["b1stab"] and LumiPageInfo[LS]["b2stab"]):
|
| 1399 |
Passed=False
|
| 1400 |
subsystemfailed.append("Beam")
|
| 1401 |
else:
|
| 1402 |
Passed=True
|
| 1403 |
|
| 1404 |
if not data_clean or (
|
| 1405 |
Rates[print_trigger]["physics"][iterator] == 1
|
| 1406 |
and Rates[print_trigger]["active"][iterator] == 1
|
| 1407 |
and Rates[print_trigger]["deadtime"][iterator] < max_dt
|
| 1408 |
#and Rates[print_trigger]["psi"][iterator] > 0
|
| 1409 |
and Passed
|
| 1410 |
and (realvalue >0.6*prediction and realvalue<1.5*prediction)
|
| 1411 |
and Rates[print_trigger]["rawrate"][iterator] > 0.04
|
| 1412 |
):
|
| 1413 |
if (print_info and num_ls==1 and (realvalue <0.4*prediction or realvalue>2.5*prediction)):
|
| 1414 |
pass
|
| 1415 |
##print '%-50s%10s%10s%10s%10s%10s%10s%10s%15s%20s' % (print_trigger,"Passed", 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)
|
| 1416 |
return True
|
| 1417 |
else:
|
| 1418 |
if (print_info and print_trigger==trig_list[0] and num_ls==1 and first_trigger):
|
| 1419 |
prediction_match=True
|
| 1420 |
if (realvalue >0.6*prediction and realvalue<1.5*prediction):
|
| 1421 |
prediction_match=False
|
| 1422 |
print '%10s%10s%10s%10s%10s%10s%10s%15s%20s%15s' % ("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, prediction_match )
|
| 1423 |
|
| 1424 |
return False
|
| 1425 |
|
| 1426 |
|
| 1427 |
#### LumiRangeGreens ####
|
| 1428 |
####inputs: RefMoreLumiArray --dict with lumi page info in LS by LS blocks,
|
| 1429 |
#### LRange --list range over lumis,
|
| 1430 |
#### nls --number of lumisections
|
| 1431 |
#### RefRunNum --run number
|
| 1432 |
####
|
| 1433 |
####outputs RangeMoreLumi --lumi page info in dict LSRange blocks with lumi, added items Run and LSRange
|
| 1434 |
def LumiRangeGreens(RefMoreLumiArray,LSRange,nls,RefRunNum,deadtimebeamactive):
|
| 1435 |
|
| 1436 |
RangeMoreLumi={}
|
| 1437 |
for keys,values in RefMoreLumiArray.iteritems():
|
| 1438 |
RangeMoreLumi[keys]=1
|
| 1439 |
|
| 1440 |
for iterator in LSRange[nls]:
|
| 1441 |
for keys, values in RefMoreLumiArray.iteritems():
|
| 1442 |
if RefMoreLumiArray[keys][iterator]==0:
|
| 1443 |
RangeMoreLumi[keys]=0
|
| 1444 |
RangeMoreLumi['LSRange']=LSRange[nls]
|
| 1445 |
RangeMoreLumi['Run']=RefRunNum
|
| 1446 |
RangeMoreLumi['DeadTimeBeamActive']=deadtimebeamactive
|
| 1447 |
return RangeMoreLumi
|
| 1448 |
|
| 1449 |
#### CheckLumis ####
|
| 1450 |
####inputs:
|
| 1451 |
#### PageLumiInfo --dict of LS with dict of some lumipage info
|
| 1452 |
#### Rates --dict of triggernames with dict of info
|
| 1453 |
def checkLS(Rates, PageLumiInfo,trig_list):
|
| 1454 |
rateslumis=Rates[trig_list[-1]]["ls"]
|
| 1455 |
keys=PageLumiInfo.keys()
|
| 1456 |
print "lumi run=",PageLumiInfo[keys[-1]]["Run"]
|
| 1457 |
ll=0
|
| 1458 |
for ls in keys:
|
| 1459 |
print ls,rateslumis[ll]
|
| 1460 |
ll=ll+1
|
| 1461 |
return False
|
| 1462 |
|
| 1463 |
|
| 1464 |
def checkL1seedChangeALLPScols(trig_list,HLTL1PS):
|
| 1465 |
|
| 1466 |
nps=0
|
| 1467 |
HLTL1_seedchanges={}
|
| 1468 |
|
| 1469 |
for HLTkey in trig_list:
|
| 1470 |
if HLTkey=='HLT_Stream_A':
|
| 1471 |
continue
|
| 1472 |
#print HLTkey
|
| 1473 |
if not HLTkey.startswith('HLT'):
|
| 1474 |
nps=9
|
| 1475 |
HLTL1_seedchanges[HLTkey]=[[0, 1, 2, 3, 4, 5, 6, 7, 8]]
|
| 1476 |
continue
|
| 1477 |
|
| 1478 |
try:
|
| 1479 |
dict=HLTL1PS[StripVersion(HLTkey)]
|
| 1480 |
#print "dict=",dict
|
| 1481 |
except:
|
| 1482 |
print "%s appears in the trigger list but does not exist in the HLT menu and is being skipped." % (StripVersion(HLTkey),)
|
| 1483 |
continue
|
| 1484 |
|
| 1485 |
HLTL1dummy={}
|
| 1486 |
for L1seed in dict.iterkeys():
|
| 1487 |
#print L1seed
|
| 1488 |
dummyL1seedlist=[]
|
| 1489 |
#print dict[L1seed]
|
| 1490 |
dummy=dict[L1seed]
|
| 1491 |
L1seedchangedummy=[]
|
| 1492 |
L1fulldummy=[]
|
| 1493 |
nps=len(dict[L1seed])
|
| 1494 |
#print "nps=",nps
|
| 1495 |
for PScol in range(0,len(dict[L1seed])):
|
| 1496 |
PScoldummy=PScol+1
|
| 1497 |
#print "PScoldummy=",PScoldummy
|
| 1498 |
if PScoldummy>(len(dict[L1seed])-1):
|
| 1499 |
PScoldummy=len(dict[L1seed])-1
|
| 1500 |
#print "changed PScoldummy=",PScoldummy
|
| 1501 |
#print PScol, PScoldummy, dummy[PScol]
|
| 1502 |
|
| 1503 |
if dummy[PScol]==dummy[PScoldummy]:
|
| 1504 |
#print PScol, "same"
|
| 1505 |
L1seedchangedummy.append(PScol)
|
| 1506 |
else:
|
| 1507 |
#print PScol, PScoldummy, "diff", dummy[PScol], dummy[PScoldummy]
|
| 1508 |
L1seedchangedummy.append(PScol)
|
| 1509 |
for ps in L1seedchangedummy:
|
| 1510 |
L1fulldummy.append(L1seedchangedummy)
|
| 1511 |
#print "L1seed change ", L1seedchangedummy, "full=",L1fulldummy
|
| 1512 |
L1seedchangedummy=[]
|
| 1513 |
for ps in L1seedchangedummy:
|
| 1514 |
L1fulldummy.append(L1seedchangedummy)
|
| 1515 |
#print "L1full=",L1fulldummy
|
| 1516 |
HLTL1dummy[L1seed]=L1fulldummy
|
| 1517 |
#print HLTL1dummy
|
| 1518 |
HLTL1_seedchanges[HLTkey]=commonL1PS(HLTL1dummy,nps)
|
| 1519 |
#print HLTkey, HLTL1_seedchanges[HLTkey]
|
| 1520 |
return HLTL1_seedchanges,nps
|
| 1521 |
|
| 1522 |
|
| 1523 |
def commonL1PS(HLTL1dummy, nps):
|
| 1524 |
### find commmon elements in L1 seeds
|
| 1525 |
HLTL1_seedchanges=[]
|
| 1526 |
for PScol in range(0,nps):
|
| 1527 |
|
| 1528 |
L1seedslist=HLTL1dummy.keys()
|
| 1529 |
L1tupletmp=set(tuple(HLTL1dummy[L1seedslist.pop()][PScol]))
|
| 1530 |
while len(L1seedslist)>0:
|
| 1531 |
L1tupletmp2=set(tuple(HLTL1dummy[L1seedslist.pop()][PScol]))
|
| 1532 |
L1tupletmp=L1tupletmp & L1tupletmp2
|
| 1533 |
if sorted(list(tuple(L1tupletmp))) not in HLTL1_seedchanges:
|
| 1534 |
HLTL1_seedchanges.append(sorted(list(tuple(L1tupletmp))))
|
| 1535 |
#print HLTL1_seedchanges
|
| 1536 |
return HLTL1_seedchanges
|
| 1537 |
|
| 1538 |
def Fitter(gr1, VX, VY, sloperate, nlow, Rates, print_trigger, first_trigger, varX, varY, lowrate):
|
| 1539 |
|
| 1540 |
f1a = 0
|
| 1541 |
f1b = 0
|
| 1542 |
f1c = 0
|
| 1543 |
f1d = 0
|
| 1544 |
if "rate" in varY:
|
| 1545 |
f1d = TF1("f1d","pol1",0,8000)#linear
|
| 1546 |
f1d.SetParameters(0.01,min(sum(VY)/sum(VX),sloperate)) ##Set Y-intercept near 0, slope either mean_rate/mean_lumi or est. slope (may be negative)
|
| 1547 |
f1d.SetLineColor(4)
|
| 1548 |
f1d.SetLineWidth(2)
|
| 1549 |
if nlow>0:
|
| 1550 |
f1d.SetParLimits(0,0,1.5*lowrate/nlow) ##Keep Y-intercept in range of low-lumi rate points
|
| 1551 |
else:
|
| 1552 |
f1d.SetParLimits(0,0,1.5*sum(VY)/len(VY))
|
| 1553 |
if (sloperate > 0):
|
| 1554 |
if (sloperate > 0.5*sum(VY)/sum(VX)): ##Slope is substantially positive
|
| 1555 |
f1d.SetParLimits(1,min(0.5*sloperate,0.5*sum(VY)/sum(VX)),1.5*sum(VY)/sum(VX))
|
| 1556 |
else: ##Slope is somewhat positive or flat
|
| 1557 |
f1d.SetParLimits(1,-0.1*sloperate,1.5*sum(VY)/sum(VX))
|
| 1558 |
else: ##Slope is negative or flat
|
| 1559 |
f1d.SetParLimits(1,1.5*sloperate,-0.1*sloperate)
|
| 1560 |
|
| 1561 |
gr1.Fit("f1d","QN","rob=0.90")
|
| 1562 |
|
| 1563 |
f1a = TF1("f1a","pol2",0,8000)#quadratic
|
| 1564 |
f1a.SetParameters(f1d.GetParameter(0),f1d.GetParameter(1),0) ##Initial values from linear fit
|
| 1565 |
f1a.SetLineColor(6)
|
| 1566 |
f1a.SetLineWidth(2)
|
| 1567 |
if nlow>0 and sloperate < 0.5*sum(VY)/sum(VX): ##Slope is not substantially positive
|
| 1568 |
f1a.SetParLimits(0,0,1.5*lowrate/nlow) ##Keep Y-intercept in range of low-lumi rate points
|
| 1569 |
else:
|
| 1570 |
f1a.SetParLimits(0,0,max(min(VY),0.3*sum(VY)/len(VY))) ##Keep Y-intercept reasonably low
|
| 1571 |
f1a.SetParLimits(1,-2.0*(max(VY)-min(VY))/(max(VX)-min(VX)),2.0*(max(VY)-min(VY))/(max(VX)-min(VX))) ##Reasonable bounds
|
| 1572 |
f1a.SetParLimits(2,-2.0*max(VY)/(max(VX)*max(VX)),2.0*max(VY)/(max(VX)*max(VX))) ##Reasonable bounds
|
| 1573 |
gr1.Fit("f1a","QN","rob=0.90")
|
| 1574 |
|
| 1575 |
if True:
|
| 1576 |
f1b = TF1("f1b","pol3",0,8000)#cubic
|
| 1577 |
f1b.SetParameters(f1a.GetParameter(0),f1a.GetParameter(1),f1a.GetParameter(2),0) ##Initial values from quadratic fit
|
| 1578 |
f1b.SetLineColor(2)
|
| 1579 |
f1b.SetLineWidth(2)
|
| 1580 |
f1b.SetParLimits(0,0,max(min(VY),0.3*sum(VY)/len(VY))) ##Keep Y-intercept reasonably low
|
| 1581 |
f1b.SetParLimits(1,-2.0*(max(VY)-min(VY))/(max(VX)-min(VX)),2.0*(max(VY)-min(VY))/(max(VX)-min(VX))) ##Reasonable bounds
|
| 1582 |
f1b.SetParLimits(2,-2.0*max(VY)/(max(VX)*max(VX)),2.0*max(VY)/(max(VX)*max(VX))) ##Reasonable bounds
|
| 1583 |
f1b.SetParLimits(3,0,2.0*max(VY)/(max(VX)*max(VX)*max(VX))) ##Reasonable bounds
|
| 1584 |
gr1.Fit("f1b","QN","rob=0.90")
|
| 1585 |
|
| 1586 |
f1c = TF1("f1c","[0]+[1]*expo(2)",0,8000)
|
| 1587 |
f1c.SetLineColor(3)
|
| 1588 |
f1c.SetLineWidth(2)
|
| 1589 |
#f1c.SetParLimits(0,0,max(min(VY),0.3*sum(VY)/len(VY)))
|
| 1590 |
f1c.SetParLimits(0,0,max(min(VY),0.01*sum(VY)/len(VY))) ##Exponential fits should start low
|
| 1591 |
f1c.SetParLimits(1,max(VY)/math.exp(15.0),max(VY)/math.exp(2.0))
|
| 1592 |
f1c.SetParLimits(2,0.0,0.0000000001)
|
| 1593 |
f1c.SetParLimits(3,2.0/max(VX),15.0/max(VX))
|
| 1594 |
gr1.Fit("f1c","QN","rob=0.90")
|
| 1595 |
##Some fits are so exponential, the graph ends early and returns a false low Chi2 value
|
| 1596 |
|
| 1597 |
else: ##If this is not a rate plot
|
| 1598 |
f1a = TF1("f1a","pol1",0,8000)
|
| 1599 |
f1a.SetLineColor(4)
|
| 1600 |
f1a.SetLineWidth(2)
|
| 1601 |
if "xsec" in varY:
|
| 1602 |
f1a.SetParLimits(0,0,meanxsec*1.5)
|
| 1603 |
if slopexsec > 0:
|
| 1604 |
f1a.SetParLimits(1,0,max(VY)/max(VX))
|
| 1605 |
else:
|
| 1606 |
f1a.SetParLimits(1,2*slopexsec,-2*slopexsec)
|
| 1607 |
else:
|
| 1608 |
f1a.SetParLimits(0,-1000,1000)
|
| 1609 |
gr1.Fit("f1a","Q","rob=0.80")
|
| 1610 |
|
| 1611 |
if (first_trigger):
|
| 1612 |
print '%-50s %4s x0 x1 x2 x3 chi2 ndf chi2/ndf' % ('trigger', 'type')
|
| 1613 |
first_trigger=False
|
| 1614 |
try:
|
| 1615 |
print '%-50s | line | % .2f | +/-%.2f | % .2e | +/-%.1e | % .2e | +/-%.1e | % .2e | +/-%.1e | %7.0f | %4.0f | %5.2f | ' % (print_trigger, f1a.GetParameter(0), f1a.GetParError(0), f1a.GetParameter(1), f1a.GetParError(1), 0 , 0 , 0 , 0 , f1a.GetChisquare(), f1a.GetNDF(), f1a_GetChisquare()/f1a.GetNDF())
|
| 1616 |
except:
|
| 1617 |
pass
|
| 1618 |
|
| 1619 |
return [f1a,f1b,f1c,f1d,first_trigger]
|
| 1620 |
|
| 1621 |
def more_fit_info(f1a,f1b,f1c,f1d,VX,VY,print_trigger,Rates):
|
| 1622 |
|
| 1623 |
meanps = median(Rates[print_trigger]["ps"])
|
| 1624 |
av_rte = mean(VY)
|
| 1625 |
passed=1
|
| 1626 |
#except ZeroDivisionError:
|
| 1627 |
try:
|
| 1628 |
f1a_Chi2 = f1a.GetChisquare()/f1a.GetNDF()
|
| 1629 |
f1b_Chi2 = f1b.GetChisquare()/f1b.GetNDF()
|
| 1630 |
f1c_Chi2 = f1c.GetChisquare()/f1c.GetNDF()
|
| 1631 |
f1d_Chi2 = f1d.GetChisquare()/f1d.GetNDF()
|
| 1632 |
except ZeroDivisionError:
|
| 1633 |
print "Zero DOF for", print_trigger
|
| 1634 |
passed=0
|
| 1635 |
f1a_BadMinimum = (f1a.GetMinimumX(5,7905,10)>2000 and f1a.GetMinimumX(5,7905,10)<7000) ##Don't allow minimum between 2000 and 7000
|
| 1636 |
f1b_BadMinimum = (f1b.GetMinimumX(5,7905,10)>2000 and f1b.GetMinimumX(5,7905,10)<7000)
|
| 1637 |
f1c_BadMinimum = ((f1c.GetMinimumX(5,7905,10)>2000 and f1c.GetMinimumX(5,7905,10)<7000)) or f1c.GetMaximum(min(VX),max(VX),10)/max(VY) > 2.0
|
| 1638 |
|
| 1639 |
return [f1a_Chi2, f1b_Chi2, f1c_Chi2,f1d_Chi2, f1a_BadMinimum, f1b_BadMinimum, f1c_BadMinimum, meanps, av_rte, passed]
|
| 1640 |
|
| 1641 |
def output_fit_info(do_fit,f1a,f1b,f1c,f1d,varX,varY,VX,VY,linear,print_trigger,first_trigger,Rates,width,chioffset,wp_bool,num_ls,meanrawrate,OutputFit, failed_paths, PSColslist, dummyPSColslist):
|
| 1642 |
[f1a_Chi2, f1b_Chi2, f1c_Chi2,f1d_Chi2, f1a_BadMinimum, f1b_BadMinimum, f1c_BadMinimum, meanps, av_rte,passed]=more_fit_info(f1a,f1b,f1c,f1d,VX,VY,print_trigger,Rates)
|
| 1643 |
OutputFit[print_trigger] = {}
|
| 1644 |
|
| 1645 |
if not do_fit:
|
| 1646 |
failure_comment= "Can't have save_fits = True and do_fit = False"
|
| 1647 |
[OutputFit,first_trigger]
|
| 1648 |
failed_paths.append([print_trigger+str(PSColslist),failure_comment])
|
| 1649 |
OutputFit[print_trigger] = ["fit failed",failure_comment]
|
| 1650 |
return [OutputFit,first_trigger]
|
| 1651 |
if min([f1a_Chi2,f1b_Chi2,f1c_Chi2,f1d_Chi2]) > 500:#require a minimum chi^2/nDOF of 500
|
| 1652 |
failure_comment = "There were events for these paths in the runs specified during the creation of the fit file, but the fit failed to converge"
|
| 1653 |
failed_paths.append([print_trigger+str(PSColslist),failure_comment])
|
| 1654 |
OutputFit[print_trigger] = ["fit failed",failure_comment]
|
| 1655 |
return [OutputFit,first_trigger]
|
| 1656 |
if "rate" in varY and not linear:
|
| 1657 |
if first_trigger:
|
| 1658 |
print '\n%-*s | TYPE | %-8s | %-11s | %-7s | %-10s | %-7s | %-10s | %-8s | %-10s | %-6s | %-4s |%-7s| %-6s |' % (width,"TRIGGER", "X0","X0 ERROR","X1","X1 ERROR","X2","X2 ERROR","X3","X3 ERROR","CHI^2","DOF","CHI2/DOF","PScols")
|
| 1659 |
first_trigger = False
|
| 1660 |
if ((f1c_Chi2 < (f1a_Chi2*chioffset) or f1a_BadMinimum) and ((f1c_Chi2 < f1b_Chi2) or f1b_BadMinimum) and f1c_Chi2 < (f1d_Chi2*chioffset) and not f1c_BadMinimum and len(VX)>1):
|
| 1661 |
graph_fit_type="expo"
|
| 1662 |
[f1c,OutputFit]=graph_output_info(f1c,graph_fit_type,print_trigger,width,num_ls,VX,VY,meanrawrate,OutputFit,PSColslist,dummyPSColslist)
|
| 1663 |
priot(wp_bool,print_trigger,meanps,f1d,f1c,graph_fit_type,av_rte)
|
| 1664 |
elif ((f1b_Chi2 < (f1a_Chi2*chioffset) or f1a_BadMinimum) and f1b_Chi2 < (f1d_Chi2*chioffset) and not f1b_BadMinimum and len(VX)>1):
|
| 1665 |
graph_fit_type="cube"
|
| 1666 |
[f1b,OutputFit]=graph_output_info(f1b,graph_fit_type,print_trigger,width,num_ls,VX,VY,meanrawrate,OutputFit,PSColslist,dummyPSColslist)
|
| 1667 |
priot(wp_bool,print_trigger,meanps,f1d,f1b,graph_fit_type,av_rte)
|
| 1668 |
elif (f1a_Chi2 < (f1d_Chi2*chioffset)):
|
| 1669 |
graph_fit_type="quad"
|
| 1670 |
[f1a,OutputFit]=graph_output_info(f1a,graph_fit_type,print_trigger,width,num_ls,VX,VY,meanrawrate,OutputFit,PSColslist,dummyPSColslist)
|
| 1671 |
priot(wp_bool,print_trigger,meanps,f1d,f1a,graph_fit_type,av_rte)
|
| 1672 |
else:
|
| 1673 |
graph_fit_type="line"
|
| 1674 |
[f1d,OutputFit]=graph_output_info(f1d,graph_fit_type,print_trigger,width,num_ls,VX,VY,meanrawrate,OutputFit,PSColslist,dummyPSColslist)
|
| 1675 |
priot(wp_bool,print_trigger,meanps,f1d,f1d,graph_fit_type,av_rte)
|
| 1676 |
elif "rate" in varY and linear:
|
| 1677 |
if first_trigger:
|
| 1678 |
print '\n%-*s | TYPE | %-8s | %-11s | %-7s | %-10s | %-7s | %-10s | %-8s | %-10s | %-6s | %-4s |%-7s| %-6s |' % (width,"TRIGGER", "X0","X0 ERROR","X1","X1 ERROR","X2","X2 ERROR","X3","X3 ERROR","CHI^2","DOF","CHI2/DOF","PScols")
|
| 1679 |
first_trigger = False
|
| 1680 |
graph_fit_type="line"
|
| 1681 |
[f1d,OutputFit]=graph_output_info(f1d,graph_fit_type,print_trigger,width,num_ls,VX,VY,meanrawrate,OutputFit,PSColslist,dummyPSColslist)
|
| 1682 |
priot(wp_bool,print_trigger,meanps,f1d,f1d,graph_fit_type,av_rte)
|
| 1683 |
else:
|
| 1684 |
graph_fit_type="quad"
|
| 1685 |
[f1a,OutputFit]=graph_output_info(f1a,graph_fit_type,print_trigger,width,num_ls,VX,VY,meanrawrate,OutputFit,PSColslist,dummyPSColslist)
|
| 1686 |
#priot(wp_bool,print_trigger,meanps,f1d,f1a,"quad",av_rte)
|
| 1687 |
|
| 1688 |
return [OutputFit,first_trigger, failed_paths]
|
| 1689 |
|
| 1690 |
def graph_output_info(graph1,graph_fit_type,print_trigger,width,num_ls,VX, VY,meanrawrate,OutputFit,PSColslist,dummyPSColslist):
|
| 1691 |
PSlist=deque(PSColslist)
|
| 1692 |
PSmin=PSlist.popleft()
|
| 1693 |
if not PSlist:
|
| 1694 |
PSmax=PSmin
|
| 1695 |
else:
|
| 1696 |
PSmax=PSlist.pop()
|
| 1697 |
|
| 1698 |
print '%-*s | %s | %-8.1f | +/-%-8.1f | %8.1e | +/-%.1e | %8.1e | +/-%.1e | %-8.1e | +/-%.1e | %6.0f | %4.0f | %5.2f | %d-%d' % (width,print_trigger, graph_fit_type,graph1.GetParameter(0) , graph1.GetParError(0) , graph1.GetParameter(1) , graph1.GetParError(1) , graph1.GetParameter(2), graph1.GetParError(2) ,graph1.GetParameter(3), graph1.GetParError(3) ,graph1.GetChisquare() , graph1.GetNDF() , graph1.GetChisquare()/graph1.GetNDF(), PSmin, PSmax)
|
| 1699 |
graph1.SetLineColor(1)
|
| 1700 |
#priot(wp_bool,print_trigger,meanps,f1d,f1c,"expo",av_rte)
|
| 1701 |
do_high_lumi = print_trigger.startswith('HLT_') and ((len(dummyPSColslist)==1 or ( max(PSColslist)>=5 and min(PSColslist)==3) ))
|
| 1702 |
sigma = CalcSigma(VX, VY, graph1, do_high_lumi)*math.sqrt(num_ls)
|
| 1703 |
OutputFit[print_trigger] = [graph_fit_type, graph1.GetParameter(0) , graph1.GetParameter(1) , graph1.GetParameter(2) ,graph1.GetParameter(3) , sigma , meanrawrate, graph1.GetParError(0) , graph1.GetParError(1) , graph1.GetParError(2) , graph1.GetParError(3)]
|
| 1704 |
|
| 1705 |
return [graph1,OutputFit]
|
| 1706 |
|
| 1707 |
def DrawFittedCurve(f1a, f1b,f1c, f1d, chioffset,do_fit,c1,VX,VY,print_trigger,Rates):
|
| 1708 |
[f1a_Chi2, f1b_Chi2, f1c_Chi2,f1d_Chi2, f1a_BadMinimum, f1b_BadMinimum, f1c_BadMinimum, meanps, av_rte, passed]=more_fit_info(f1a,f1b,f1c,f1d,VX,VY,print_trigger,Rates)
|
| 1709 |
|
| 1710 |
|
| 1711 |
if do_fit:
|
| 1712 |
try:
|
| 1713 |
if ((f1c_Chi2 < (f1a_Chi2*chioffset) or f1a_BadMinimum ) and (f1c_Chi2 < f1b_Chi2 or f1b_BadMinimum ) and not f1c_BadMinimum ):
|
| 1714 |
f1c.Draw("same")
|
| 1715 |
elif ( (f1b_Chi2 < (f1a_Chi2*chioffset) or f1a_BadMinimum) and not f1b_BadMinimum):
|
| 1716 |
f1b.Draw("same")
|
| 1717 |
else:
|
| 1718 |
f1a.Draw("same")
|
| 1719 |
|
| 1720 |
f1d.Draw("same")
|
| 1721 |
except:
|
| 1722 |
True
|
| 1723 |
|
| 1724 |
c1.Update()
|
| 1725 |
|
| 1726 |
return c1
|
| 1727 |
|
| 1728 |
def EndMkrootfile(failed_paths, save_fits, save_root, fit_file, RootFile, OutputFit,OutputFitPS,L1SeedChangeFit):
|
| 1729 |
|
| 1730 |
if len(failed_paths) > 0:
|
| 1731 |
if save_fits:
|
| 1732 |
print "\n***************NO FIT RECORDED FOR THE FOLLOWING PATHS***************"
|
| 1733 |
else:
|
| 1734 |
print "\n***************THE FOLLOWING PATHS HAVE BEEN SKIPPED BECAUSE THE FIT WAS MISSING***************"
|
| 1735 |
sorted_failed_paths = sorted(failed_paths, key=itemgetter(1))
|
| 1736 |
for error_comment, entries in groupby(sorted_failed_paths, key=itemgetter(1)):
|
| 1737 |
print '\n'+error_comment+':'
|
| 1738 |
if 'not enough datapoints' in error_comment:
|
| 1739 |
print "(For a given trigger, if a group of PS columns has been skipped, the fit to all PS columns will be used in that region.)"
|
| 1740 |
for entry in entries:
|
| 1741 |
print entry[0]
|
| 1742 |
|
| 1743 |
if save_root:
|
| 1744 |
print "\nOutput root file is "+str(RootFile)
|
| 1745 |
#print "DONE:",OutputFit
|
| 1746 |
if save_fits:
|
| 1747 |
if os.path.exists(fit_file):
|
| 1748 |
os.remove(fit_file)
|
| 1749 |
FitOutputFile = open(fit_file, 'wb')
|
| 1750 |
pickle.dump(OutputFit, FitOutputFile, 2)
|
| 1751 |
FitOutputFile.close()
|
| 1752 |
print "Output fit file is "+str(fit_file)
|
| 1753 |
if save_fits and L1SeedChangeFit:
|
| 1754 |
PSfitfile=fit_file.replace("HLT_NoV","HLT_NoV_ByPS")
|
| 1755 |
print "A corresponding PS fit file has been saved."
|
| 1756 |
if os.path.exists(PSfitfile):
|
| 1757 |
os.remove(PSfitfile)
|
| 1758 |
FitOutputFilePS= open(PSfitfile, 'wb')
|
| 1759 |
pickle.dump(OutputFitPS,FitOutputFilePS,2)
|
| 1760 |
FitOutputFilePS.close()
|
| 1761 |
|
| 1762 |
##### NEED BETTER gr1 def for failure#####
|
| 1763 |
def DefineGraphs(print_trigger,OutputFit,do_fit,varX,varY,x_label,y_label,VX,VY,VXE,VYE,VF,VFE,fit_file, failed_paths,PSColslist):
|
| 1764 |
passed=1
|
| 1765 |
try:
|
| 1766 |
gr1 = TGraphErrors(len(VX), VX, VY, VXE, VYE)
|
| 1767 |
|
| 1768 |
except:
|
| 1769 |
failure_comment = "In runs specified during creation of the fit file, there were no events for this path: probably due to high deadtime or low raw (prescaled) rate"
|
| 1770 |
failed_paths.append([print_trigger,failure_comment])
|
| 1771 |
if do_fit:
|
| 1772 |
OutputFit[print_trigger] = ["fit failed",failure_comment]
|
| 1773 |
#gr1 = TGraphErrors(1, VX, VY, VXE, VYE)
|
| 1774 |
#gr3 = TGraphErrors(1, VX, VF, VXE, VFE)
|
| 1775 |
###replaces continue in main fucntion
|
| 1776 |
passed=0
|
| 1777 |
return [OutputFit,0, 0, failed_paths, passed]
|
| 1778 |
try:
|
| 1779 |
if not do_fit:
|
| 1780 |
gr3 = TGraphErrors(len(VX), VX, VF, VXE, VFE)
|
| 1781 |
else:
|
| 1782 |
##fake defn (will not be used)
|
| 1783 |
gr3 =TGraphErrors(len(VX), VX, VY, VXE, VYE)
|
| 1784 |
except:
|
| 1785 |
print "gr3 failed to define!"
|
| 1786 |
|
| 1787 |
exit(2)
|
| 1788 |
|
| 1789 |
|
| 1790 |
if not do_fit:
|
| 1791 |
gr3.SetMarkerStyle(8)
|
| 1792 |
gr3.SetMarkerSize(0.4)
|
| 1793 |
gr3.SetMarkerColor(4)
|
| 1794 |
gr3.SetFillColor(4)
|
| 1795 |
gr3.SetFillStyle(3003)
|
| 1796 |
|
| 1797 |
|
| 1798 |
if (len(VX)<10 and do_fit):
|
| 1799 |
failure_comment = "In runs specified during creation of the fit file, there were not enough datapoints for these paths (probably due to high deadtime or low raw (prescaled) rate)"
|
| 1800 |
failed_paths.append([print_trigger+" PS columns: "+str(PSColslist),failure_comment])
|
| 1801 |
OutputFit[print_trigger] = ["fit failed",failure_comment]
|
| 1802 |
gr1 = TGraphErrors(1, VX, VY, VXE, VYE)
|
| 1803 |
###replaces continue in main fucntion
|
| 1804 |
passed=0
|
| 1805 |
return [OutputFit,gr1, gr3,failed_paths, passed]
|
| 1806 |
|
| 1807 |
gr1.SetName("Graph_"+str(print_trigger)+"_"+str(varY)+"_vs_"+str(varX))
|
| 1808 |
gr1.GetXaxis().SetTitle(x_label)
|
| 1809 |
gr1.GetYaxis().SetTitle(y_label)
|
| 1810 |
gr1.SetTitle(str(print_trigger))
|
| 1811 |
gr1.SetMinimum(0)
|
| 1812 |
gr1.SetMaximum(1.2*max(VY))
|
| 1813 |
#gr1.GetXaxis().SetLimits(min(VX)-0.2*max(VX),1.2*max(VX))
|
| 1814 |
gr1.GetXaxis().SetLimits(0,1.2*max(VX))
|
| 1815 |
gr1.SetMarkerStyle(8)
|
| 1816 |
|
| 1817 |
if fit_file:
|
| 1818 |
gr1.SetMarkerSize(0.8)
|
| 1819 |
else:
|
| 1820 |
gr1.SetMarkerSize(0.5)
|
| 1821 |
gr1.SetMarkerColor(2)
|
| 1822 |
|
| 1823 |
|
| 1824 |
return [OutputFit,gr1, gr3, failed_paths, passed]
|
| 1825 |
|
| 1826 |
def DrawSave(save_root, save_png, var, varY, print_trigger, do_fit, gr1, gr3, chioffset, f1a, f1b, f1c, f1d, RootFile):
|
| 1827 |
if save_root or save_png:
|
| 1828 |
c1 = TCanvas(str(varX),str(varY))
|
| 1829 |
c1.SetName(str(print_trigger)+"_"+str(varY)+"_vs_"+str(varX))
|
| 1830 |
gr1.Draw("APZ")
|
| 1831 |
if not do_fit:
|
| 1832 |
gr3.Draw("P3")
|
| 1833 |
c1.Update()
|
| 1834 |
else:
|
| 1835 |
c1=DrawFittedCurve(f1a, f1b, f1c, f1d, chioffset,do_fit,c1)
|
| 1836 |
|
| 1837 |
if save_root:
|
| 1838 |
myfile = TFile( RootFile, 'UPDATE' )
|
| 1839 |
c1.Write()
|
| 1840 |
myfile.Close()
|
| 1841 |
if save_png:
|
| 1842 |
c1.SaveAs(str(print_trigger)+"_"+str(varY)+"_vs_"+str(varX)+".png")
|
| 1843 |
|
| 1844 |
|
| 1845 |
if __name__=='__main__':
|
| 1846 |
global thisyear
|
| 1847 |
main()
|
