1313
1414#Make sure process is single-threaded, we'll use multiprocessing library instead
1515os .environ ["OMP_NUM_THREADS" ] = "1"
16+ nCores = multiprocessing .cpu_count ()- 1
1617
1718# Reduces warning messages.
1819ROOT .RooMsgService .instance ().setGlobalKillBelow (ROOT .RooFit .WARNING )
154155expectedBackgroundCounts = []
155156for i in range (0 ,nSources ):
156157 expectedBackgroundCounts .append (bgndObject .countsInFitRange * sourceObjects [i ].normalization / float (bgndObject .normalization ))
158+ #TODO: dead time calculator.
159+ #source run rate - bgnd run rate = increased source rate.
160+ #Increased source rate * waveform length = fractional dead time
161+ #Expected bgnd rates (w/dead time) = expectedBackgroundCounts*(1-fractional dead time)
162+
157163expectedSourceCounts = []
158164for i in range (0 ,nSources ):
159165 expectedSourceCounts .append (sourceObjects [i ].countsInFitRange - expectedBackgroundCounts [i ])
@@ -240,9 +246,11 @@ def lnlike(theta):
240246 ROOT .RooFit .Verbose (0 ),
241247 ROOT .RooFit .NumCPU (1 )
242248 )
249+ ROOT .SetOwnership (nll ,True )
243250
244251 #Make NLL positive
245252 nllVal += (- 1 * nll .getVal ())
253+
246254
247255 if PLOT == 1 :
248256 modelHist = hist .Clone ("modelHist" )
@@ -254,11 +262,7 @@ def lnlike(theta):
254262 bgndHist = ROOT .RooAbsData .createHistogram (bgndDataHist ,"energyVar" ,energyVar ,ROOT .RooFit .Binning ("binning" ))
255263 bgndCounts = expectedBackgroundCounts [i ]
256264
257- print (expectedSourceCounts [i ])
258- print (sourceScaling )
259265 simCounts = sourceScaling * expectedSourceCounts [i ]
260- print (sourceObjects [i ].name ,sourceCounts ,bgndCounts ,simCounts )
261-
262266 helperFns .plotResults (sourceHist ,sourceCounts ,bgndHist ,bgndCounts ,tempHist ,simCounts ,modelHist ,title + "_fit_" + sourceObjects [i ].name )
263267 del sourceDataHist
264268 del sourceHist
@@ -290,7 +294,7 @@ def lnprob(theta):
290294
291295
292296##########################MC RUN STARTS HERE####################################
293- with multiprocessing .Pool () as pool :
297+ with multiprocessing .Pool (nCores ) as pool :
294298 sampler = emcee .EnsembleSampler (nWalkers ,ndim ,lnprob ,pool = pool )
295299 #Burn in
296300 print ("Starting burn in..." )
@@ -308,7 +312,7 @@ def lnprob(theta):
308312samples = sampler .flatchain
309313lnprobs = sampler .lnprobability [:,:]
310314flatLnprobs = lnprobs .reshape (- 1 )
311- with open (title + "_sampler .csv" , 'w' ) as samplerOutFile :
315+ with open ("sampler_" + title + ".csv" , 'w' ) as samplerOutFile :
312316 theWriter = csvlib .writer (samplerOutFile , delimiter = ',' )
313317 for sampleLine , llval in zip (samples , flatLnprobs ):
314318 theWriter .writerow (np .append (sampleLine ,llval ))
@@ -321,12 +325,12 @@ def lnprob(theta):
321325 axes = fig .add_subplot (gs [i ,:])
322326 axes .plot (sampler .chain [:,:,i ].T , '-' , color = 'k' , alpha = 0.3 )
323327 axes .set_title (labels [i ])
324- plt .savefig (title + "_traces .pdf" )
328+ plt .savefig ("tracePlot_" + title + ".pdf" )
325329
326330#CORNER PLOT HERE
327331samples = sampler .flatchain
328332fig = corner .corner (samples , labels = labels , ranges = ranges , quantiles = [0.16 ,0.5 ,0.84 ],show_titles = True ,title_kwargs = {'fontsize' :12 })
329- fig .savefig (title + "corner .pdf" )
333+ fig .savefig ("corner_" + title + ".pdf" )
330334
331335#CALCULATE QUANTILES HERE
332336bestFitValues = []
0 commit comments