Thoughts

Should IAV effect and NL (quenching/Cerenkov) be folded together? Will this improve the uncertainty on the lowest bins? Technically the combined effect should vary depending on whether the positron passes through the IAV at the beginning or at the end of its range.

What’s a super histogram?

From comment above FluxCalculator::LoadSuperHistograms:

Super-histograms are the neutrino flux (binned in true neutrino energy) seen in one AD in an entire period (.e.g. 20 weeks), using the weekly efficiencies and livetimes to weigh each element in the sum.

Questions for Henoch

[X] What’s up with “x+2” in ReactorTheoreticalPrediction_P15A_8AD.py? Related to taking the midpoint of the 5 bins we avg over
[ ] Why do we only use the “D1” (?) SNF spectrum?
[X] Should we change line 176 of make_combined_spectra.C to:

	    + snf_fraction * flux_snf[0][i] * get_ibd_xsec(enu[i]) /flux_snf_total[0]*flux_wo_snf_total_ave;

Actually Yasu’s version seems correct since the above will give us the number of observed (i.e. xsec-weighted) events. We report the raw flux, so need to divide out the xsec.

[ ] Difference between reactor/fissionIsotopeSpectra_Huber_v0.txt and DBS/p15a_reactor/Files/EnChuan_fissionIsotopeSpectra_Huber_Fit.txt Former not used?

Generating tags

ctags-exuberant -e -R .

Misc

OscProbTable (and hence Ranger) are only(?) referenced by ShapeFit/fit_shape_sterile_hybrid_scan.C

As well as by ShapeFit/LoadClasses.C

How is PredictedIBD generated?

Inputs:

ShapeFit/Flux/superhistograms : gen by genSuperHistograms.C, used by FluxCalculator ShapeFit/PredictedIBD/PredictedIBD_asdf.root : gen by genToySpectraTree???? used by Predictor

Reactor weekly information

See Henoch’s mywork/DayaBaySpectrum/p15a_reactor (WeeklyAvg.dat, Matt’s weekly_livetime files) see also p14a_uinblind8AD_isotope – Yasu files The output files from ^ go into Theta13Analysis/reactor_covmatrix (Unblinded6AD.txt etc.) Which is used by quickFit/reactor/Spectrum.C

Where is SNF and NE added in? reactor_covmatrix/make_combined_spectrum.C (has Wei’s numbers from DayaBaySpectrum/reactor) cf. reactor_covmatrix/p15a/isotope_spectra

What we need to add to Doxygen?

reactor_covmatrix? DayaBaySpectrum/p15a_reactor?

My current understanding

Raw inputs:

weekly_livetime (from Matt) –> DayaBaySpectrum/p15a_reactor –> reactor_covmatrix/Unblinded6AD.txt etc. –>

Henoch’s email

Update ToySpectra as you want
Generate toy MC in order to generate covariance matrix
1. open genToySpectraTree.C and check the background files are correct
2. open rungenToySpectraTree.C and double check the input and

destination files

execute “root rungenToySpectraTree.C”

Generate evis-to-enu matrix

***Note: So far I’ve been using the evis-to-enu matrix generated from the P12E period***

matrix_evis_to_enu_rateonly_unified_fine_p12e_unblinded_noosc.txt (Rate only fit) matrix_evis_to_enu_unified_fine_p12e_unblinded.txt (Rate + Shape fit) matrix_evis_to_enu_unified_fine_p12e_unblinded_noosc.txt (For generating Covariance matrix)

(Not needed unless flux input file from Christine is changed)
1. execute “genSuperHistograms.C”
Generate covariance matrices
1. Double check background file and the default oscillation

parameters in build_covmatrix.C

Double check input and output file names in run_build_covmatrix.C
Execute run_build_covmatrix.sh
cd covariance_matrices

By default this should create a covariance matrix for all signal systematic, and a covariance matrix for all background systematic.

Fit!
1. The current fit code is:

Rate+Shape: fit_shape_2d_6.C fit_shape_2d_8.C fit_shape_2d_6_8.C

Rate Only: fit_rateonly_6.C fit_rateonly_8.C fit_rateonly_6_8.C

Double check all input files…. .
1. Super histogram
2. Theta13 input file (txt file)
3. Background file
4. Enu to evis conversion file
run whatever fit macros….

The 6, 8 and 6_8 corresponds doing the fit on the 6AD period only, 8AD period only and 6AD+8AD periods combined. The codes are very similar, in calculating the chi2 there is a flag where 0=6AD, 1=8AD, -1=6AD+8AD.

Log of commands for P15A test-fit

Prepare flux predictions

Uncorrected prediction

In DayaBaySpectrum/p15a_reactor:

python ReactorTheoreticalPrediction_P15A_6AD.py
python ReactorTheoreticalPrediction_P15A_8AD.py

Output (in output/; symlinked to /reactor_covmatrix/p15a/isotope_spectra):

reactor_P15A_6AD.txt reactor_P15A_6AD_U235.txt reactor_P15A_8AD_Pu239.txt reactor_P15A_8AD_U238.txt reactor_P15A_6AD_Pu239.txt reactor_P15A_6AD_U238.txt reactor_P15A_8AD_Pu241.txt reactor_P15A_6AD_Pu241.txt reactor_P15A_8AD.txt reactor_P15A_8AD_U235.txt

Prediction w/ NE/SNF

In reactor_covmatrix/p15a:

root -b -q make_combined_spectra.C

Output (here):

reactor_P15A_6AD_SNF_nonEq.txt reactor_P15A_8AD_SNF_nonEq.txt

Generate “data files” to configure toy MC

Already done, but if needed: In toySpectra/data_file:

./generate_data_file.py
# ./generate_data_file_extra.py

Note: The “extra” configs, which don’t seem to be in use, are reactor_only, det_only, allsys_wo_reactor, allsys_wo_det, and {nominal,allsys,sigsys,reactor_only,allsys_wo_det}_bcwflux

Generate toy samples

We need the sigsys and bgsys samples in order to generate the covmatrix. I had to edit runGenToySpectraTree to enable bgsys In toySpectra:

root -l 'rungenToySpectraTree.C(3)'
root -l 'rungenToySpectraTree.C(4)'

Note: Used my “old” ROOT 5.34.38. CINT segfaults at the end because yeah.

Generate evis-to-enu matrix

Seems like we just need matrix_evis_to_enu_fine.txt. Had to modify genEvisToEnuMatrix.C to restore LBNL binning from “MaximRequest” binning. In toySpectra, within ROOT (export LBNL_FIT_DEBUG=1):

root -b -q LoadClasses.C genEvisToEnuMatrix.C+

Then, in ShapeFit:

root -b -q make_evis_to_enu_matrix_fine.C

Generate super histograms

In toySpectra:

root -b -q LoadClasses.C genSuperHistograms.C+

Generate PredictedIBD

In toySpectra:

root -b -q LoadClasses.C genPredictedIBD.C+

Build covariance matrices

In ShapeFit:

root -b -q 'run_build_covmatrix.C(9)' # sigsys
root -b -q 'run_build_covmatrix.C(18)' # bgsys

NOTE There will be scary “Ay caramba” warnings about the discrepancies between the spectra rates (from toy MC) and text rates (from data). This is fine; it’s just the reactor antineutrino anomaly (???). If these warnings appear when we’re actually performing a fit, then we have a real issue.

Fit!

In ShapeFit:

root -b -q LoadClasses.C fit_shape_2d_P15A.C+

Stuff I had to copy over

…from Henoch to the “pristine” tree (henochw.bak)

For Spectrum.C

Defined in Spectrum.C

reactor_covmatrix/p15a/nNu_Mcov_combined_huber-french_u238cor.txt
toySpectra/unified_nl_data/consModel_450itr.root

Defined in data_file

toySpectra/IavDistortion/iavMatrix_P14A.root

In genToySpectraTree.C

toySpectra/unblinded_baseline.txt
ShapeFit/Inputs/Theta13-inputs_P15A_inclusive_6ad.txt
ShapeFit/Inputs/Theta13-inputs_P15A_inclusive_8ad.txt
ShapeFit/Spectra/accidental_eprompt_shapes_6ad.root
ShapeFit/Spectra/accidental_eprompt_shapes_8ad.root

For fit_shape_2d_P15A.C

ShapeFit/Spectra/ibd_eprompt_shapes_6ad_LBNL.root
ShapeFit/Spectra/ibd_eprompt_shapes_8ad_LBNL.root

Stuff to provide for a new fit:

ShapeFit/Inputs/Theta13-inputs_P17B_inclusive_8ad.txt
ShapeFit/Spectra/accidental_eprompt_shapes_8ad.root
ShapeFit/Spectra/ibd_eprompt_shapes_8ad_LBNL.root
toySpectra/data_file/dyb_data_v1_nominal.txt (if generating a new covmatrix at updated theta/delta)

Things to check

Differences in osc params for noosc vs nominal data files
Why in genToySpectraTree, we multiply by binWidth for signal spectrum but not for bkg? Looks OK, see test_binning.C. Bkg spectrum is a true histogram whereas Sig spectrum uses /MeV units. Should confirm in Spectrum.C.

Log of commands for P17B test-fit

Done in “beda.test”, cloned from “beda.mine”.

Prepare directories

scripts/prep_dirs.sh

Provide input files

“Theta13” file (sig/bkg rates, efficiencies, livetimes, …)
IBD spectra (see Beda’s prep scripts in ShapeFit/Spectra)
Accidental spectra (ditto)
Replace dyb_data_v1_nominal.txt (if generating a new covmatrix)

P17B example inputs

scripts/install_example.sh

Set environment variables

For official (IHEP) result, need to enable BCW binning:

export LBNL_FIT_BINNING=BCW

Prepare flux predictions

XXX Update day-array sizes in Produce_Isotope_Spectra_XYZ.C

Uncorrected prediction

In ReactorPowerCalculator: XXX We don’t know whether to use BCW core spectra (seemingly the default for Beda) or the LBNL spectra (as for P15A). I personally (re)set the default to LBNL for now; the below will use LBNL. Update: Toy MC doesn’t like the LBNL binning, using BCW option=1 below.

root -b -q 'Produce_Isotope_SpectraP17B_unblinded.C(1)'

Prediction w/ NE/SNF

In ReactorPowerCalculator/isotope_spectra_by_Beda:

root -b -q make_combined_spectra_P17B_unblinded.C

Generate “data files” to configure toy MC

In toySpectra/data_file:

./generate_data_file.py
# ./generate_data_file_extra.py

Note: The “extra” configs, which don’t seem to be in use, are reactor_only, det_only, allsys_wo_reactor, allsys_wo_det, and {nominal,allsys,sigsys,reactor_only,allsys_wo_det}_bcwflux

Generate toy samples

We need the sigsys and bgsys samples in order to generate the covmatrix. I uncommented allsys just to get the numbers to line up. In toySpectra:

root -l 'rungenToySpectraTree.C(2)'
root -l 'rungenToySpectraTree.C(3)'

Generate evis-to-enu matrix

Seems like we just need matrix_evis_to_enu_fine.txt. Had to modify genEvisToEnuMatrix.C to restore LBNL binning (37 evis) from BCW binning (26 evis). In toySpectra:

root -b -q LoadClasses.C genEvisToEnuMatrix.C+

Then, in ShapeFit:

root -b -q make_evis_to_enu_matrix_fine_P17B.C

Generate super histograms

In toySpectra:

root -b -q LoadClasses.C genSuperHistograms.C+

Generate PredictedIBD

In toySpectra:

root -b -q LoadClasses.C genPredictedIBD.C+

XXX Should we use rungenPredictedIBD instead?

Build covariance matrices

In ShapeFit:

root -b -q 'run_build_covmatrix.C(9)' # sigsys
root -b -q 'run_build_covmatrix.C(21)' # bgsys

NOTE There will be scary “Ay caramba” warnings about the discrepancies between the spectra rates (from toy MC) and text rates (from data). This is fine; it’s just the reactor antineutrino anomaly (???). If these warnings appear when we’re actually performing a fit, then we have a real issue.

Fit!

In ShapeFit:

root -b -q LoadClasses.C fit_shape_2d_P17B.C+

Fiterating

Set t13/dm2 in data_file
Set t13 for genPredictedIBD
…?

Things to consider

Consistent random seed?

Headscratchers

Why does genPredictedIBD use a theta of zero whereas genToySpectraTree uses nominal?
What were the consequences of genPredictedIBD using fine evis binning previously?

Benchmarking test_chain.sh

Whole things takes 72 minutes on Cori login.
48 minutes with bg/sys parallelization

One-time steps

genReactor: 20sec
genToyConf: 1sec

Each-time steps, nToys-dependent => ~1hr @ 1000

nToys = 10

genToys: 45sec
genCovMat: 40sec

nToys = 30

genToys: 1min
genCovMat: 1min

nToys = 100

genToys 1m40s => 16min @ 1000
genCovMat: 3m30s => 35 min @ 1000

Each-time steps, nToys-independent => 15min

genEvisEnu: 15min
genSuperHists: 20sec
genPredIBD: 15sec

The fit, grid-size (θ x Δ) dependent

Default is 31 x 31

(10 x 10)

shapeFit: 3min

(20 x 20)

shapeFit: 8min

(31 x 31)

shapeFit: 18min

Things to parallelize

genToys (5 min; 2.5 min when sig/bg par): toys -> 1min (2procs x 10 = 20 threads, ~2000% CPU) XXX reduce number of threads? Need to re-tune?
genCovMat (30 min; 15 min when sig/bg par): toys -> 3.5 min (2procs x 8 = 16 threads, ~1200% CPU) :tofinetune:
genEvisEnu (7 min): bins -> 1.5 min (30 threads) (could keep going; currently at 21% mem)
shapeFit (17 min): gridpoints -> 8 min (12 threads) :(

Bottlenecks?

Reading reactor covariance matrix

New things to check

make_combined_spectra_P17B_new_SNF.C? See snf_coarse_percentage
Does MINUIT fail to converge when we use the LBNL or IHEP inputs? (It “fails” with our input although the result is consistent with LBNL’s)

Notes on multithreading

ROOT Forum says (next link) a TFile can only be accessed from the thread that constructed it. “Accessing” a TFile includes, of course, writing to it, even if it done indirectly, via Write methods of other objects. and Using the same TFile from different threads, even protecting concurrent access with a mutex, as well as creating a TFile in one thread and using it from another, can result in quirky behavior. Yet we do it anyway (with a mutex). Should each thread create its own TFile?
OpenMP - ROOT - ROOT Forum
Multicore/multithreading - ROOT - ROOT Forum (N. Tagg)

Calling functions from <omp.h>

If code is compiled by ACLiC, do, e.g.:

#ifdef __CLING__
extern int omp_get_thread_num();
#else
#include <omp.h>
#endif

OpenMP references

http://on-demand.gputechconf.com/gtc/2016/presentation/s6510-jeff-larkin-targeting-gpus-openmp.pdf
https://www.exascaleproject.org/wp-content/uploads/2017/05/OpenMP-4.5-and-Beyond-SOLLVE-part-21.pdf

Files

NOTES.org

Latest commit

History

NOTES.org

File metadata and controls

Thoughts

What’s a super histogram?

Questions for Henoch

Generating tags

Misc

Reactor weekly information

What we need to add to Doxygen?

My current understanding

Henoch’s email

Log of commands for P15A test-fit

Prepare flux predictions

Uncorrected prediction

Prediction w/ NE/SNF

Generate “data files” to configure toy MC

Generate toy samples

Generate evis-to-enu matrix

Generate super histograms

Generate PredictedIBD

Build covariance matrices

Fit!

Stuff I had to copy over

For Spectrum.C

Defined in Spectrum.C

Defined in data_file

In genToySpectraTree.C

For fit_shape_2d_P15A.C

Stuff to provide for a new fit:

Things to check

Log of commands for P17B test-fit

Prepare directories

Provide input files

P17B example inputs

Set environment variables

Prepare flux predictions

Uncorrected prediction

Prediction w/ NE/SNF

Generate “data files” to configure toy MC

Generate toy samples

Generate evis-to-enu matrix

Generate super histograms

Generate PredictedIBD

Build covariance matrices

Fit!

Fiterating

Things to consider

Headscratchers

Benchmarking test_chain.sh

One-time steps

Each-time steps, nToys-dependent => ~1hr @ 1000

nToys = 10

nToys = 30

nToys = 100

Each-time steps, nToys-independent => 15min

The fit, grid-size (θ x Δ) dependent

(10 x 10)

(20 x 20)

(31 x 31)

Things to parallelize

Bottlenecks?

New things to check

Notes on multithreading

Calling functions from <omp.h>

OpenMP references