diff --git a/data/MicroBooNE/CC1Mu1p/CC1Mu1p_BinScheme.txt b/data/MicroBooNE/CC1Mu1p/CC1Mu1p_BinScheme.txt new file mode 100644 index 000000000..6e14fda16 --- /dev/null +++ b/data/MicroBooNE/CC1Mu1p/CC1Mu1p_BinScheme.txt @@ -0,0 +1,426 @@ +SerialDeltaPT_DeltaAlphaTPlot + +Bin # & Slice & Low edge & High edge + +1 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0 & 0.05 +2 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.05 & 0.1 +3 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.1 & 0.15 +4 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.15 & 0.2 +5 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.2 & 0.25 +6 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.25 & 0.3 +7 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.3 & 0.35 +8 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.35 & 0.4 +9 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.4 & 0.47 +10 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.47 & 0.55 +11 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.55 & 0.9 +12 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0 & 0.05 +13 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.05 & 0.1 +14 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.1 & 0.15 +15 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.15 & 0.2 +16 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.2 & 0.25 +17 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.25 & 0.3 +18 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.3 & 0.35 +19 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.35 & 0.4 +20 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.4 & 0.47 +21 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.47 & 0.55 +22 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.55 & 0.65 +23 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.65 & 0.9 +24 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0 & 0.05 +25 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.05 & 0.1 +26 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.1 & 0.15 +27 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.15 & 0.2 +28 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.2 & 0.25 +29 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.25 & 0.3 +30 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.3 & 0.35 +31 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.35 & 0.4 +32 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.4 & 0.47 +33 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.47 & 0.55 +34 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.55 & 0.65 +35 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.65 & 0.75 +36 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.75 & 0.9 +37 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0 & 0.05 +38 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.05 & 0.1 +39 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.1 & 0.15 +40 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.15 & 0.2 +41 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.2 & 0.25 +42 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.25 & 0.3 +43 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.3 & 0.35 +44 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.35 & 0.4 +45 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.4 & 0.47 +46 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.47 & 0.55 +47 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.55 & 0.65 +48 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.65 & 0.75 +49 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.75 & 0.9 + + +SerialDeltaPT_MuonCosThetaPlot + +Bin # & Slice & Low edge & High edge + +1 & -1 < cos#theta_{#mu} < 0 & 0 & 0.05 +2 & -1 < cos#theta_{#mu} < 0 & 0.05 & 0.1 +3 & -1 < cos#theta_{#mu} < 0 & 0.1 & 0.15 +4 & -1 < cos#theta_{#mu} < 0 & 0.15 & 0.2 +5 & -1 < cos#theta_{#mu} < 0 & 0.2 & 0.25 +6 & -1 < cos#theta_{#mu} < 0 & 0.25 & 0.3 +7 & -1 < cos#theta_{#mu} < 0 & 0.3 & 0.35 +8 & -1 < cos#theta_{#mu} < 0 & 0.35 & 0.4 +9 & -1 < cos#theta_{#mu} < 0 & 0.4 & 0.47 +10 & -1 < cos#theta_{#mu} < 0 & 0.47 & 0.55 +11 & -1 < cos#theta_{#mu} < 0 & 0.55 & 0.65 +12 & -1 < cos#theta_{#mu} < 0 & 0.65 & 0.75 +13 & -1 < cos#theta_{#mu} < 0 & 0.75 & 0.9 +14 & 0 < cos#theta_{#mu} < 0.5 & 0 & 0.05 +15 & 0 < cos#theta_{#mu} < 0.5 & 0.05 & 0.1 +16 & 0 < cos#theta_{#mu} < 0.5 & 0.1 & 0.15 +17 & 0 < cos#theta_{#mu} < 0.5 & 0.15 & 0.2 +18 & 0 < cos#theta_{#mu} < 0.5 & 0.2 & 0.25 +19 & 0 < cos#theta_{#mu} < 0.5 & 0.25 & 0.3 +20 & 0 < cos#theta_{#mu} < 0.5 & 0.3 & 0.35 +21 & 0 < cos#theta_{#mu} < 0.5 & 0.35 & 0.4 +22 & 0 < cos#theta_{#mu} < 0.5 & 0.4 & 0.47 +23 & 0 < cos#theta_{#mu} < 0.5 & 0.47 & 0.55 +24 & 0 < cos#theta_{#mu} < 0.5 & 0.55 & 0.65 +25 & 0 < cos#theta_{#mu} < 0.5 & 0.65 & 0.75 +26 & 0 < cos#theta_{#mu} < 0.5 & 0.75 & 0.9 +27 & 0.5 < cos#theta_{#mu} < 0.75 & 0 & 0.05 +28 & 0.5 < cos#theta_{#mu} < 0.75 & 0.05 & 0.1 +29 & 0.5 < cos#theta_{#mu} < 0.75 & 0.1 & 0.15 +30 & 0.5 < cos#theta_{#mu} < 0.75 & 0.15 & 0.2 +31 & 0.5 < cos#theta_{#mu} < 0.75 & 0.2 & 0.25 +32 & 0.5 < cos#theta_{#mu} < 0.75 & 0.25 & 0.3 +33 & 0.5 < cos#theta_{#mu} < 0.75 & 0.3 & 0.35 +34 & 0.5 < cos#theta_{#mu} < 0.75 & 0.35 & 0.4 +35 & 0.5 < cos#theta_{#mu} < 0.75 & 0.4 & 0.47 +36 & 0.5 < cos#theta_{#mu} < 0.75 & 0.47 & 0.55 +37 & 0.5 < cos#theta_{#mu} < 0.75 & 0.55 & 0.65 +38 & 0.5 < cos#theta_{#mu} < 0.75 & 0.65 & 0.75 +39 & 0.5 < cos#theta_{#mu} < 0.75 & 0.75 & 0.9 +40 & 0.75 < cos#theta_{#mu} < 1 & 0 & 0.05 +41 & 0.75 < cos#theta_{#mu} < 1 & 0.05 & 0.1 +42 & 0.75 < cos#theta_{#mu} < 1 & 0.1 & 0.15 +43 & 0.75 < cos#theta_{#mu} < 1 & 0.15 & 0.2 +44 & 0.75 < cos#theta_{#mu} < 1 & 0.2 & 0.25 +45 & 0.75 < cos#theta_{#mu} < 1 & 0.25 & 0.3 +46 & 0.75 < cos#theta_{#mu} < 1 & 0.3 & 0.35 +47 & 0.75 < cos#theta_{#mu} < 1 & 0.35 & 0.4 +48 & 0.75 < cos#theta_{#mu} < 1 & 0.4 & 0.47 +49 & 0.75 < cos#theta_{#mu} < 1 & 0.47 & 0.55 +50 & 0.75 < cos#theta_{#mu} < 1 & 0.55 & 0.65 +51 & 0.75 < cos#theta_{#mu} < 1 & 0.65 & 0.75 +52 & 0.75 < cos#theta_{#mu} < 1 & 0.75 & 0.9 + + +SerialDeltaPT_ProtonCosThetaPlot + +Bin # & Slice & Low edge & High edge + +1 & -1 < cos#theta_{p} < 0 & 0 & 0.1 +2 & -1 < cos#theta_{p} < 0 & 0.1 & 0.2 +3 & -1 < cos#theta_{p} < 0 & 0.2 & 0.3 +4 & -1 < cos#theta_{p} < 0 & 0.3 & 0.4 +5 & -1 < cos#theta_{p} < 0 & 0.4 & 0.55 +6 & -1 < cos#theta_{p} < 0 & 0.55 & 0.65 +7 & -1 < cos#theta_{p} < 0 & 0.65 & 0.75 +8 & -1 < cos#theta_{p} < 0 & 0.75 & 0.9 +9 & 0 < cos#theta_{p} < 0.5 & 0 & 0.05 +10 & 0 < cos#theta_{p} < 0.5 & 0.05 & 0.1 +11 & 0 < cos#theta_{p} < 0.5 & 0.1 & 0.15 +12 & 0 < cos#theta_{p} < 0.5 & 0.15 & 0.2 +13 & 0 < cos#theta_{p} < 0.5 & 0.2 & 0.25 +14 & 0 < cos#theta_{p} < 0.5 & 0.25 & 0.3 +15 & 0 < cos#theta_{p} < 0.5 & 0.3 & 0.35 +16 & 0 < cos#theta_{p} < 0.5 & 0.35 & 0.4 +17 & 0 < cos#theta_{p} < 0.5 & 0.4 & 0.47 +18 & 0 < cos#theta_{p} < 0.5 & 0.47 & 0.55 +19 & 0 < cos#theta_{p} < 0.5 & 0.55 & 0.65 +20 & 0 < cos#theta_{p} < 0.5 & 0.65 & 0.75 +21 & 0 < cos#theta_{p} < 0.5 & 0.75 & 0.9 +22 & 0.5 < cos#theta_{p} < 0.75 & 0 & 0.05 +23 & 0.5 < cos#theta_{p} < 0.75 & 0.05 & 0.1 +24 & 0.5 < cos#theta_{p} < 0.75 & 0.1 & 0.15 +25 & 0.5 < cos#theta_{p} < 0.75 & 0.15 & 0.2 +26 & 0.5 < cos#theta_{p} < 0.75 & 0.2 & 0.25 +27 & 0.5 < cos#theta_{p} < 0.75 & 0.25 & 0.3 +28 & 0.5 < cos#theta_{p} < 0.75 & 0.3 & 0.35 +29 & 0.5 < cos#theta_{p} < 0.75 & 0.35 & 0.4 +30 & 0.5 < cos#theta_{p} < 0.75 & 0.4 & 0.47 +31 & 0.5 < cos#theta_{p} < 0.75 & 0.47 & 0.55 +32 & 0.5 < cos#theta_{p} < 0.75 & 0.55 & 0.65 +33 & 0.5 < cos#theta_{p} < 0.75 & 0.65 & 0.75 +34 & 0.5 < cos#theta_{p} < 0.75 & 0.75 & 0.9 +35 & 0.75 < cos#theta_{p} < 1 & 0 & 0.05 +36 & 0.75 < cos#theta_{p} < 1 & 0.05 & 0.1 +37 & 0.75 < cos#theta_{p} < 1 & 0.1 & 0.15 +38 & 0.75 < cos#theta_{p} < 1 & 0.15 & 0.2 +39 & 0.75 < cos#theta_{p} < 1 & 0.2 & 0.25 +40 & 0.75 < cos#theta_{p} < 1 & 0.25 & 0.3 +41 & 0.75 < cos#theta_{p} < 1 & 0.3 & 0.35 +42 & 0.75 < cos#theta_{p} < 1 & 0.35 & 0.4 +43 & 0.75 < cos#theta_{p} < 1 & 0.4 & 0.47 +44 & 0.75 < cos#theta_{p} < 1 & 0.47 & 0.55 +45 & 0.75 < cos#theta_{p} < 1 & 0.55 & 0.65 +46 & 0.75 < cos#theta_{p} < 1 & 0.65 & 0.9 + + +SerialDeltaAlphaT_DeltaPTPlot + +Bin # & Slice & Low edge & High edge + +1 & 0.0 < #deltap_{T} < 0.2 GeV/c & 0 & 22 +2 & 0.0 < #deltap_{T} < 0.2 GeV/c & 22 & 44 +3 & 0.0 < #deltap_{T} < 0.2 GeV/c & 44 & 66 +4 & 0.0 < #deltap_{T} < 0.2 GeV/c & 66 & 88 +5 & 0.0 < #deltap_{T} < 0.2 GeV/c & 88 & 110 +6 & 0.0 < #deltap_{T} < 0.2 GeV/c & 110 & 145 +7 & 0.0 < #deltap_{T} < 0.2 GeV/c & 145 & 180 +8 & 0.2 < #deltap_{T} < 0.4 GeV/c & 0 & 22 +9 & 0.2 < #deltap_{T} < 0.4 GeV/c & 22 & 44 +10 & 0.2 < #deltap_{T} < 0.4 GeV/c & 44 & 66 +11 & 0.2 < #deltap_{T} < 0.4 GeV/c & 66 & 88 +12 & 0.2 < #deltap_{T} < 0.4 GeV/c & 88 & 110 +13 & 0.2 < #deltap_{T} < 0.4 GeV/c & 110 & 145 +14 & 0.2 < #deltap_{T} < 0.4 GeV/c & 145 & 180 +15 & 0.4 < #deltap_{T} < 1.0 GeV/c & 0 & 22 +16 & 0.4 < #deltap_{T} < 1.0 GeV/c & 22 & 44 +17 & 0.4 < #deltap_{T} < 1.0 GeV/c & 44 & 66 +18 & 0.4 < #deltap_{T} < 1.0 GeV/c & 66 & 88 +19 & 0.4 < #deltap_{T} < 1.0 GeV/c & 88 & 110 +20 & 0.4 < #deltap_{T} < 1.0 GeV/c & 110 & 145 +21 & 0.4 < #deltap_{T} < 1.0 GeV/c & 145 & 180 + + +SerialDeltaAlphaT_MuonCosThetaPlot + +Bin # & Slice & Low edge & High edge + +1 & -1 < cos#theta_{#mu} < 0 & 0 & 22 +2 & -1 < cos#theta_{#mu} < 0 & 22 & 44 +3 & -1 < cos#theta_{#mu} < 0 & 44 & 66 +4 & -1 < cos#theta_{#mu} < 0 & 66 & 88 +5 & -1 < cos#theta_{#mu} < 0 & 88 & 110 +6 & -1 < cos#theta_{#mu} < 0 & 110 & 145 +7 & -1 < cos#theta_{#mu} < 0 & 145 & 180 +8 & 0 < cos#theta_{#mu} < 0.5 & 0 & 22 +9 & 0 < cos#theta_{#mu} < 0.5 & 22 & 44 +10 & 0 < cos#theta_{#mu} < 0.5 & 44 & 66 +11 & 0 < cos#theta_{#mu} < 0.5 & 66 & 88 +12 & 0 < cos#theta_{#mu} < 0.5 & 88 & 110 +13 & 0 < cos#theta_{#mu} < 0.5 & 110 & 145 +14 & 0 < cos#theta_{#mu} < 0.5 & 145 & 180 +15 & 0.5 < cos#theta_{#mu} < 0.75 & 0 & 22 +16 & 0.5 < cos#theta_{#mu} < 0.75 & 22 & 44 +17 & 0.5 < cos#theta_{#mu} < 0.75 & 44 & 66 +18 & 0.5 < cos#theta_{#mu} < 0.75 & 66 & 88 +19 & 0.5 < cos#theta_{#mu} < 0.75 & 88 & 110 +20 & 0.5 < cos#theta_{#mu} < 0.75 & 110 & 145 +21 & 0.5 < cos#theta_{#mu} < 0.75 & 145 & 180 +22 & 0.75 < cos#theta_{#mu} < 1 & 0 & 22 +23 & 0.75 < cos#theta_{#mu} < 1 & 22 & 44 +24 & 0.75 < cos#theta_{#mu} < 1 & 44 & 66 +25 & 0.75 < cos#theta_{#mu} < 1 & 66 & 88 +26 & 0.75 < cos#theta_{#mu} < 1 & 88 & 110 +27 & 0.75 < cos#theta_{#mu} < 1 & 110 & 145 +28 & 0.75 < cos#theta_{#mu} < 1 & 145 & 180 + + +SerialDeltaAlphaT_ProtonCosThetaPlot + +Bin # & Slice & Low edge & High edge + +1 & -1 < cos#theta_{p} < 0 & 0 & 22 +2 & -1 < cos#theta_{p} < 0 & 22 & 44 +3 & -1 < cos#theta_{p} < 0 & 44 & 66 +4 & -1 < cos#theta_{p} < 0 & 66 & 88 +5 & -1 < cos#theta_{p} < 0 & 88 & 110 +6 & -1 < cos#theta_{p} < 0 & 110 & 145 +7 & -1 < cos#theta_{p} < 0 & 145 & 180 +8 & 0 < cos#theta_{p} < 0.5 & 0 & 22 +9 & 0 < cos#theta_{p} < 0.5 & 22 & 44 +10 & 0 < cos#theta_{p} < 0.5 & 44 & 66 +11 & 0 < cos#theta_{p} < 0.5 & 66 & 88 +12 & 0 < cos#theta_{p} < 0.5 & 88 & 110 +13 & 0 < cos#theta_{p} < 0.5 & 110 & 145 +14 & 0 < cos#theta_{p} < 0.5 & 145 & 180 +15 & 0.5 < cos#theta_{p} < 0.75 & 0 & 22 +16 & 0.5 < cos#theta_{p} < 0.75 & 22 & 44 +17 & 0.5 < cos#theta_{p} < 0.75 & 44 & 66 +18 & 0.5 < cos#theta_{p} < 0.75 & 66 & 88 +19 & 0.5 < cos#theta_{p} < 0.75 & 88 & 110 +20 & 0.5 < cos#theta_{p} < 0.75 & 110 & 145 +21 & 0.5 < cos#theta_{p} < 0.75 & 145 & 180 +22 & 0.75 < cos#theta_{p} < 1 & 0 & 22 +23 & 0.75 < cos#theta_{p} < 1 & 22 & 44 +24 & 0.75 < cos#theta_{p} < 1 & 44 & 66 +25 & 0.75 < cos#theta_{p} < 1 & 66 & 88 +26 & 0.75 < cos#theta_{p} < 1 & 88 & 110 +27 & 0.75 < cos#theta_{p} < 1 & 110 & 145 +28 & 0.75 < cos#theta_{p} < 1 & 145 & 180 + + +SerialDeltaPhiT_DeltaPTPlot + +Bin # & Slice & Low edge & High edge + +1 & 0.0 < #deltap_{T} < 0.2 GeV/c & 0 & 12.5 +2 & 0.0 < #deltap_{T} < 0.2 GeV/c & 12.5 & 25 +3 & 0.0 < #deltap_{T} < 0.2 GeV/c & 25 & 37.5 +4 & 0.0 < #deltap_{T} < 0.2 GeV/c & 37.5 & 180 +5 & 0.2 < #deltap_{T} < 0.4 GeV/c & 0 & 20 +6 & 0.2 < #deltap_{T} < 0.4 GeV/c & 20 & 40 +7 & 0.2 < #deltap_{T} < 0.4 GeV/c & 40 & 60 +8 & 0.2 < #deltap_{T} < 0.4 GeV/c & 60 & 80 +9 & 0.2 < #deltap_{T} < 0.4 GeV/c & 80 & 100 +10 & 0.2 < #deltap_{T} < 0.4 GeV/c & 100 & 120 +11 & 0.2 < #deltap_{T} < 0.4 GeV/c & 120 & 140 +12 & 0.2 < #deltap_{T} < 0.4 GeV/c & 140 & 160 +13 & 0.2 < #deltap_{T} < 0.4 GeV/c & 160 & 180 +14 & 0.4 < #deltap_{T} < 1.0 GeV/c & 0 & 25 +15 & 0.4 < #deltap_{T} < 1.0 GeV/c & 25 & 50 +16 & 0.4 < #deltap_{T} < 1.0 GeV/c & 50 & 75 +17 & 0.4 < #deltap_{T} < 1.0 GeV/c & 75 & 105 +18 & 0.4 < #deltap_{T} < 1.0 GeV/c & 105 & 145 +19 & 0.4 < #deltap_{T} < 1.0 GeV/c & 145 & 180 + + +SerialDeltaPtx_DeltaPtyPlot + +Bin # & Slice & Low edge & High edge + +1 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & -0.55 & -0.45 +2 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & -0.45 & -0.35 +3 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & -0.35 & -0.25 +4 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & -0.25 & -0.15 +5 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & -0.15 & -0.05 +6 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & -0.05 & 0.05 +7 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.05 & 0.15 +8 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.15 & 0.25 +9 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.25 & 0.35 +10 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.35 & 0.45 +11 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.45 & 0.55 +12 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & -0.55 & -0.45 +13 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & -0.45 & -0.35 +14 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & -0.35 & -0.25 +15 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & -0.25 & -0.15 +16 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & -0.15 & -0.05 +17 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & -0.05 & 0.05 +18 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.05 & 0.15 +19 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.15 & 0.25 +20 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.25 & 0.35 +21 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.35 & 0.45 +22 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.45 & 0.55 +23 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & -0.55 & -0.35 +24 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & -0.35 & -0.25 +25 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & -0.25 & -0.15 +26 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & -0.15 & -0.05 +27 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & -0.05 & 0.05 +28 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.05 & 0.15 +29 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.15 & 0.25 +30 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.25 & 0.35 +31 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.35 & 0.55 + + +SerialECal_DeltaPTPlot + +Bin # & Slice & Low edge & High edge + +1 & 0.0 < #deltap_{T} < 0.2 GeV/c & 0.2 & 0.34 +2 & 0.0 < #deltap_{T} < 0.2 GeV/c & 0.34 & 0.48 +3 & 0.0 < #deltap_{T} < 0.2 GeV/c & 0.48 & 0.62 +4 & 0.0 < #deltap_{T} < 0.2 GeV/c & 0.62 & 0.76 +5 & 0.0 < #deltap_{T} < 0.2 GeV/c & 0.76 & 0.9 +6 & 0.0 < #deltap_{T} < 0.2 GeV/c & 0.9 & 1.04 +7 & 0.0 < #deltap_{T} < 0.2 GeV/c & 1.04 & 1.18 +8 & 0.0 < #deltap_{T} < 0.2 GeV/c & 1.18 & 1.32 +9 & 0.0 < #deltap_{T} < 0.2 GeV/c & 1.32 & 1.6 +10 & 0.2 < #deltap_{T} < 0.4 GeV/c & 0.2 & 0.34 +11 & 0.2 < #deltap_{T} < 0.4 GeV/c & 0.34 & 0.48 +12 & 0.2 < #deltap_{T} < 0.4 GeV/c & 0.48 & 0.62 +13 & 0.2 < #deltap_{T} < 0.4 GeV/c & 0.62 & 0.76 +14 & 0.2 < #deltap_{T} < 0.4 GeV/c & 0.76 & 0.9 +15 & 0.2 < #deltap_{T} < 0.4 GeV/c & 0.9 & 1.04 +16 & 0.2 < #deltap_{T} < 0.4 GeV/c & 1.04 & 1.18 +17 & 0.2 < #deltap_{T} < 0.4 GeV/c & 1.18 & 1.32 +18 & 0.2 < #deltap_{T} < 0.4 GeV/c & 1.32 & 1.6 +19 & 0.4 < #deltap_{T} < 1.0 GeV/c & 0.2 & 0.34 +20 & 0.4 < #deltap_{T} < 1.0 GeV/c & 0.34 & 0.48 +21 & 0.4 < #deltap_{T} < 1.0 GeV/c & 0.48 & 0.62 +22 & 0.4 < #deltap_{T} < 1.0 GeV/c & 0.62 & 0.76 +23 & 0.4 < #deltap_{T} < 1.0 GeV/c & 0.76 & 0.9 +24 & 0.4 < #deltap_{T} < 1.0 GeV/c & 0.9 & 1.04 +25 & 0.4 < #deltap_{T} < 1.0 GeV/c & 1.04 & 1.18 +26 & 0.4 < #deltap_{T} < 1.0 GeV/c & 1.18 & 1.32 +27 & 0.4 < #deltap_{T} < 1.0 GeV/c & 1.32 & 1.6 + + +SerialECal_DeltaAlphaTPlot + +Bin # & Slice & Low edge & High edge + +1 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.2 & 0.38 +2 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.38 & 0.55 +3 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.55 & 0.73 +4 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.73 & 0.9 +5 & 0^{o} < #delta#alpha_{T} < 45^{o} & 0.9 & 1.08 +6 & 0^{o} < #delta#alpha_{T} < 45^{o} & 1.08 & 1.23 +7 & 0^{o} < #delta#alpha_{T} < 45^{o} & 1.23 & 1.4 +8 & 0^{o} < #delta#alpha_{T} < 45^{o} & 1.4 & 1.6 +9 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.2 & 0.34 +10 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.34 & 0.48 +11 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.48 & 0.62 +12 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.62 & 0.76 +13 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.76 & 0.9 +14 & 45^{o} < #delta#alpha_{T} < 90^{o} & 0.9 & 1.04 +15 & 45^{o} < #delta#alpha_{T} < 90^{o} & 1.04 & 1.18 +16 & 45^{o} < #delta#alpha_{T} < 90^{o} & 1.18 & 1.39 +17 & 45^{o} < #delta#alpha_{T} < 90^{o} & 1.39 & 1.6 +18 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.2 & 0.34 +19 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.34 & 0.48 +20 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.48 & 0.62 +21 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.62 & 0.76 +22 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.76 & 0.9 +23 & 90^{o} < #delta#alpha_{T} < 135^{o} & 0.9 & 1.04 +24 & 90^{o} < #delta#alpha_{T} < 135^{o} & 1.04 & 1.18 +25 & 90^{o} < #delta#alpha_{T} < 135^{o} & 1.18 & 1.39 +26 & 90^{o} < #delta#alpha_{T} < 135^{o} & 1.39 & 1.6 +27 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.2 & 0.38 +28 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.38 & 0.55 +29 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.55 & 0.73 +30 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.73 & 0.9 +31 & 135^{o} < #delta#alpha_{T} < 180^{o} & 0.9 & 1.08 +32 & 135^{o} < #delta#alpha_{T} < 180^{o} & 1.08 & 1.25 +33 & 135^{o} < #delta#alpha_{T} < 180^{o} & 1.25 & 1.43 +34 & 135^{o} < #delta#alpha_{T} < 180^{o} & 1.43 & 1.6 + + +SerialECal_DeltaPtyPlot + +Bin # & Slice & Low edge & High edge + +1 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.2 & 0.34 +2 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.34 & 0.48 +3 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.48 & 0.62 +4 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.62 & 0.76 +5 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.76 & 0.9 +6 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 0.9 & 1.04 +7 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 1.04 & 1.18 +8 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 1.18 & 1.32 +9 & -0.75 < #deltap_{T,y} < -0.15 GeV/c & 1.32 & 1.6 +10 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.2 & 0.34 +11 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.34 & 0.48 +12 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.48 & 0.62 +13 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.62 & 0.76 +14 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.76 & 0.9 +15 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 0.9 & 1.04 +16 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 1.04 & 1.18 +17 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 1.18 & 1.32 +18 & -0.15 < #deltap_{T,y} < 0.15 GeV/c & 1.32 & 1.6 +19 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.2 & 0.34 +20 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.34 & 0.48 +21 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.48 & 0.62 +22 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.62 & 0.76 +23 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.76 & 0.9 +24 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 0.9 & 1.04 +25 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 1.04 & 1.18 +26 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 1.18 & 1.32 +27 & 0.15 < #deltap_{T,y} < 0.45 GeV/c & 1.32 & 1.6 diff --git a/data/MicroBooNE/CC1Mu1p/CC1Mu1p_DataRelease.root b/data/MicroBooNE/CC1Mu1p/CC1Mu1p_DataRelease.root new file mode 100644 index 000000000..4951751a3 Binary files /dev/null and b/data/MicroBooNE/CC1Mu1p/CC1Mu1p_DataRelease.root differ diff --git a/src/FCN/SampleList.cxx b/src/FCN/SampleList.cxx index 0eab78a45..7ec48c3cb 100644 --- a/src/FCN/SampleList.cxx +++ b/src/FCN/SampleList.cxx @@ -186,6 +186,7 @@ #include "MicroBooNE_CC1MuNp_XSec_1D_nu.h" #include "MicroBooNE_CC1ENp_XSec_1D_nu.h" #include "MicroBooNE_CCInc_XSec_2DPcos_nu.h" +#include "MicroBooNE_CC1Mu1p_XSec_2D_nu.h" #endif #ifdef MINERvA_ENABLED @@ -1146,7 +1147,19 @@ MeasurementBase *CreateSample(nuiskey samplekey) { !name.compare("MicroBooNE_CC1Mu1p_XSec_1DECal_nu") || !name.compare("MicroBooNE_CC1Mu1p_XSec_1DEQE_nu")) { return (new MicroBooNE_CC1Mu1p_XSec_1D_nu(samplekey)); - } else + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPT_DeltaAlphaT_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPT_MuonCosTheta_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPT_ProtonCosTheta_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaAlphaT_DeltaPT_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaAlphaT_MuonCosTheta_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaAlphaT_ProtonCosTheta_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPhiT_DeltaPT_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPtx_DeltaPty_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DECal_DeltaPT_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DECal_DeltaAlphaT_nu") || + !name.compare("MicroBooNE_CC1Mu1p_XSec_2DECal_DeltaPty_nu")) { + return (new MicroBooNE_CC1Mu1p_XSec_2D_nu(samplekey)); + } else #endif #ifdef MINERvA_ENABLED diff --git a/src/MicroBooNE/CMakeLists.txt b/src/MicroBooNE/CMakeLists.txt index 71bc689f8..365021ba7 100644 --- a/src/MicroBooNE/CMakeLists.txt +++ b/src/MicroBooNE/CMakeLists.txt @@ -17,6 +17,7 @@ # along with NUISANCE. If not, see . ################################################################################ set(MicroBooNE_Impl_Files + MicroBooNE_CC1Mu1p_XSec_2D_nu.cxx MicroBooNE_CCInc_XSec_2DPcos_nu.cxx MicroBooNE_CC1MuNp_XSec_1D_nu.cxx MicroBooNE_CC1ENp_XSec_1D_nu.cxx @@ -26,6 +27,7 @@ set(MicroBooNE_Impl_Files ) set(MicroBooNE_Hdr_Files + MicroBooNE_CC1Mu1p_XSec_2D_nu.h MicroBooNE_CCInc_XSec_2DPcos_nu.h MicroBooNE_CC1MuNp_XSec_1D_nu.h MicroBooNE_CC1Mu2p_XSec_1D_nu.h diff --git a/src/MicroBooNE/MicroBooNE_CC1Mu1p_XSec_2D_nu.cxx b/src/MicroBooNE/MicroBooNE_CC1Mu1p_XSec_2D_nu.cxx new file mode 100644 index 000000000..580138a7b --- /dev/null +++ b/src/MicroBooNE/MicroBooNE_CC1Mu1p_XSec_2D_nu.cxx @@ -0,0 +1,702 @@ +// Copyright 2016 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret + +/******************************************************************************* + * This file is part of NUISANCE. + * + * NUISANCE is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * NUISANCE is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with NUISANCE. If not, see . + *******************************************************************************/ + +#include "MicroBooNE_CC1Mu1p_XSec_2D_nu.h" +#include "MicroBooNE_SignalDef.h" + +#include "TH1D.h" +#include "TH2D.h" + +namespace { + +// Extract numbers from string (ignore units and text) +std::vector extractNumbers(const std::string &text) { + std::vector nums; + std::regex numRegex(R"([-+]?\d*\.?\d+)"); + for (std::sregex_iterator it(text.begin(), text.end(), numRegex), end; it != end; ++it) { + nums.push_back(std::stod(it->str())); + } + return nums; +} + +// Parse the slice name string +std::pair parseSliceEdges(const std::string &slice, + const double default_low=0.0, const double default_high=1.0) { + double low = -std::numeric_limits::infinity(); + double high = std::numeric_limits::infinity(); + + // Identify patterns of comparisons + if (slice.find('<') != std::string::npos || slice.find('>') != std::string::npos) { + std::vector nums = extractNumbers(slice); + + if (nums.size() == 2 && slice.find('<') != std::string::npos) { + low = std::min(nums[0], nums[1]); + high = std::max(nums[0], nums[1]); + } else if (nums.size() == 1) { + if (slice.find('>') != std::string::npos) + low = nums[0]; + else if (slice.find('<') != std::string::npos) + high = nums[0]; + } + } + + // Apply defaults if missing + if (std::isinf(low)) low = default_low; + if (std::isinf(high)) high = default_high; + + return {low, high}; +} + +void DivideBinWidth(TH1D* h) { + + int NBins = h->GetXaxis()->GetNbins(); + + for (int i = 0; i < NBins; i++) { + + double CurrentEntry = h->GetBinContent(i+1); + double NewEntry = CurrentEntry / h->GetBinWidth(i+1); + + double CurrentError = h->GetBinError(i+1); + double NewError = CurrentError / h->GetBinWidth(i+1); + + h->SetBinContent(i+1,NewEntry); + h->SetBinError(i+1,NewError); + + } + +} + +} // namespace + +Slice::Slice(std::string name, int id) : slice_name(name), slice_id(id) { + slice_edges = parseSliceEdges(slice_name); +} + +Slice::~Slice() {} + +void Slice::AddBin(int global, double low, double high) { + // TODO: error handling + global_bins.push_back(global); + bin_edges.emplace_back(low, high); +} + +double* Slice::GetBinsForTH1() const { + int n_bins = global_bins.size(); + double *bin_edges_arr = new double[n_bins + 1]; // there's always 1 more edge than bins + for (int i = 0; i < n_bins; ++i) { + bin_edges_arr[i] = bin_edges[i].first; + } + bin_edges_arr[n_bins] = bin_edges[n_bins-1].second; + return bin_edges_arr; +} + +BinScheme::BinScheme() {} + +BinScheme::~BinScheme() {} + +void BinScheme::AddBin(std::vector bin_config) { + // Cast values as parsed by LoadBinScheme + std::string slice_name = bin_config[1]; + int global_bin = std::stoi(bin_config[0]); + double low_edge = std::stod(bin_config[2]); + double high_edge = std::stod(bin_config[3]); + return AddBin(slice_name, global_bin, low_edge, high_edge); +} + +void BinScheme::AddBin(std::string slice_name, + int global, double low, double high) { + for (auto& slice : slice_vec) { + // Add bin to slice if already exists in collection + if (slice.GetSliceName() == slice_name) { + slice.AddBin(global, low, high); + return; + } + } + // Create new slice if not found before + Slice new_slice(slice_name, current_slice); + current_slice++; + new_slice.AddBin(global, low, high); + slice_vec.push_back(new_slice); +} + +int BinScheme::GetNumberBins() { + int nbins = 0; + for (auto& slice : slice_vec) { + nbins += slice.GetNumberBins(); + } + return nbins; +} + +Slice BinScheme::GetSliceFromGlobal(int global) { + for (auto& slice : slice_vec) { + std::vector global_in_slice = slice.GetGlobalBins(); + if (std::find(global_in_slice.begin(), global_in_slice.end(), global) != global_in_slice.end()) { + return slice; + } + } + NUIS_ABORT("MicroBooNE_CC1Mu1p_XSec_2D_nu: Invalid bin number in scheme!"); +} + +//******************************************************************** +MicroBooNE_CC1Mu1p_XSec_2D_nu::MicroBooNE_CC1Mu1p_XSec_2D_nu( + nuiskey samplekey) { +//******************************************************************** + fSettings = LoadSampleSettings(samplekey); + + std::string name = fSettings.GetS("name"); + + // work out which sample you need, and set axii + if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPT_DeltaAlphaT_nu")) { + fDist = kDeltaPT; + fSlice = kDeltaAlphaT; + fSuffix = "SerialDeltaPT_DeltaAlphaT"; + fSliceTitle = "#delta#alpha_{T} (deg)"; + fSettings.SetXTitle("#deltap_{T} (GeV/c)"); + fSettings.SetYTitle("d^{2}#sigma/d#delta#alpha_{T}d#deltap_{T} (cm^{2}/(deg)/(GeV/c)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPT_MuonCosTheta_nu")) { + fDist = kDeltaPT; + fSlice = kMuonCosTheta; + fSuffix = "SerialDeltaPT_MuonCosTheta"; + fSliceTitle = "cos#theta_{#mu}"; + fSettings.SetXTitle("#deltap_{T} (GeV/c)"); + fSettings.SetYTitle("d^{2}#sigma/dcos#theta_{#mu}d#deltap_{T} (cm^{2}/(GeV/c)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPT_ProtonCosTheta_nu")) { + fDist = kDeltaPT; + fSlice = kProtonCosTheta; + fSuffix = "SerialDeltaPT_ProtonCosTheta"; + fSliceTitle = "cos#theta_{p}"; + fSettings.SetXTitle("#deltap_{T} (GeV/c)"); + fSettings.SetYTitle("d^{2}#sigma/dcos#theta_{p}d#deltap_{T} (cm^{2}/(GeV/c)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaAlphaT_DeltaPT_nu")) { + fDist = kDeltaAlphaT; + fSlice = kDeltaPT; + fSuffix = "SerialDeltaAlphaT_DeltaPT"; + fSliceTitle = "#deltap_{T} (GeV/c)"; + fSettings.SetXTitle("#delta#alpha_{T} (deg)"); + fSettings.SetYTitle("d^{2}#sigma/d#deltap_{T}d#delta#alpha_{T} (cm^{2}/(GeV/c)/(deg)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaAlphaT_MuonCosTheta_nu")) { + fDist = kDeltaAlphaT; + fSlice = kMuonCosTheta; + fSuffix = "SerialDeltaAlphaT_MuonCosTheta"; + fSliceTitle = "cos#theta_{#mu}"; + fSettings.SetXTitle("#delta#alpha_{T} (deg)"); + fSettings.SetYTitle("d^{2}#sigma/dcos#theta_{#mu}d#delta#alpha_{T} (cm^{2}/(deg)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaAlphaT_ProtonCosTheta_nu")) { + fDist = kDeltaAlphaT; + fSlice = kProtonCosTheta; + fSuffix = "SerialDeltaAlphaT_ProtonCosTheta"; + fSliceTitle = "cos#theta_{p}"; + fSettings.SetXTitle("#delta#alpha_{T} (deg)"); + fSettings.SetYTitle("d^{2}#sigma/dcos#theta_{p}d#delta#alpha_{T} (cm^{2}/(deg)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPhiT_DeltaPT_nu")) { + fDist = kDeltaPhiT; + fSlice = kDeltaPT; + fSuffix = "SerialDeltaPhiT_DeltaPT"; + fSliceTitle = "#deltap_{T} (GeV/c)"; + fSettings.SetXTitle("#delta#phi_{T} (deg)"); + fSettings.SetYTitle("d^{2}#sigma/dcos#theta_{p}d#delta#phi_{T} (cm^{2}/(deg)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DDeltaPtx_DeltaPty_nu")) { + fDist = kDeltaPtx; + fSlice = kDeltaPty; + fSuffix = "SerialDeltaPtx_DeltaPty"; + fSliceTitle = "#deltap_{T,y} (GeV/c)"; + fSettings.SetXTitle("#deltap_{T,x} (GeV/c)"); + fSettings.SetYTitle("d^{2}#sigma/d#deltap_{T,y}d#deltap_{T,x} (cm^{2}/(GeV/c)^{2}/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DECal_DeltaPT_nu")) { + fDist = kECal; + fSlice = kDeltaPT; + fSuffix = "SerialECal_DeltaPT"; + fSliceTitle = "#deltap_{T} (GeV/c)"; + fSettings.SetXTitle("E^{cal} (GeV)"); + fSettings.SetYTitle("d^{2}#sigma/d#deltap_{T}dE^{ecal} (cm^{2}/(GeV/c)/(GeV)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DECal_DeltaAlphaT_nu")) { + fDist = kECal; + fSlice = kDeltaAlphaT; + fSuffix = "SerialECal_DeltaAlphaT"; + fSliceTitle = "#delta#alpha_{T} (deg)"; + fSettings.SetXTitle("E^{cal} (GeV)"); + fSettings.SetYTitle("d^{2}#sigma/d#delta#alpha_{T}dE^{ecal} (cm^{2}/(deg)/(GeV)/^{40}Ar)"); + } else if (!name.compare("MicroBooNE_CC1Mu1p_XSec_2DECal_DeltaPty_nu")) { + fDist = kECal; + fSlice = kDeltaPty; + fSuffix = "SerialECal_DeltaPty"; + fSliceTitle = "#deltap_{T,y} (GeV/c)"; + fSettings.SetXTitle("E^{cal} (GeV)"); + fSettings.SetYTitle("d^{2}#sigma/d#deltap_{T,y}dE^{ecal} (cm^{2}/(GeV/c)/(GeV)/^{40}Ar)"); + } else { + NUIS_ABORT( + "MicroBooNE_CC1Mu1p_XSec_2D_nu: Didn't get a valid name: " << name); + } + + // Sample overview --------------------------------------------------- + std::string descrip = name + " sample.\n" + "Target: Ar\n" + "Flux: BNB FHC numu\n" + "Signal: CC1Mu1p\n"; + fSettings.SetDescription(descrip); + fSettings.SetTitle(name); + fSettings.SetAllowedTypes("FULL,DIAG/FREE,SHAPE,FIX/SYSTCOV/STATCOV", + "FIX/FULL"); + fSettings.SetEnuRange(0.0, 6.8); + fSettings.DefineAllowedTargets("Ar"); + fSettings.DefineAllowedSpecies("numu"); + FinaliseSampleSettings(); + + // Load data --------------------------------------------------------- + SetHistograms(); + + // ScaleFactor for DiffXSec/cm2/Nucleus + // Already multiplied by the Ar mass number + fScaleFactor = GetEventHistogram()->Integral("width") / fNEvents * 1.0E-38 / + TotalIntegratedFlux() * 40; + + // Final setup ------------------------------------------------------ + FinaliseMeasurement(); +} + +bool MicroBooNE_CC1Mu1p_XSec_2D_nu::isSignal(FitEvent *event) { + return SignalDef::MicroBooNE::isCC1Mu1p(event, EnuMin, EnuMax); +} + +void MicroBooNE_CC1Mu1p_XSec_2D_nu::FillEventVariables(FitEvent *event) { + + if (!isSignal(event)) { // double the work, but it lets us use the below + // functions without error checking + fXVar = -999; + return; + } + + auto const &signal_proton = + *SignalDef::MicroBooNE::GetCC1Mu1pProtonsInPS(event).front(); + TVector3 vpmu = event->GetHMFSParticle(13)->P3(); + + // using definitions in + // https://journals.aps.org/prd/pdf/10.1103/PhysRevD.108.053002 + + TVector3 vp = signal_proton.P3(); + TVector3 vSum = + vpmu + + vp; // missing momentum in the plane transverse to the beam direction + + TVector3 vpmuT(vpmu.X(), vpmu.Y(), 0.); // transverse muon momentum + TVector3 vpT(vp.X(), vp.Y(), 0.); // transverse proton momentum + TVector3 vSumT(vSum.X(), vSum.Y(), 0.); // transverse missing momentum (eq. 1) + + TVector3 vpmuL(0., 0., vpmu.Z()); + TVector3 vpL(0., 0., vp.Z()); + + double DeltaPT = vSum.Pt() / 1000.; // GeV/c + double DeltaAlphaT = + TMath::ACos((-vpmuT * vSumT) / (vpmuT.Mag() * vSumT.Mag())) * 180. / + TMath::Pi(); // deg + double DeltaPhiT = + TMath::ACos((-vpmuT * vpT) / (vpmuT.Mag() * vpT.Mag())) * 180. / + TMath::Pi(); // deg -////////////changed by Abi 5/12/23 using (eq 3) + + double MuonCosTheta = vpmu.CosTheta(); + double ProtonCosTheta = vp.CosTheta(); + + double MuonMomentum = vpmu.Mag() / 1000.; // GeV/c + double ProtonMomentum = vp.Mag() / 1000.; // GeV/c + + double BE = 0.04; // GeV, binding energy + double NeutronMass_GeV = 0.939565; // GeV + + double ProtonMass_GeV = 0.938272; // GeV + double MuonMass_GeV = 0.106; // GeV + double DeltaM2 = TMath::Power(NeutronMass_GeV, 2.) - + TMath::Power(ProtonMass_GeV, 2.); // GeV^2 + double MuonEnergy = (event->GetHMFSParticle(13)->fP.E()) / 1000.0; // GeV + // Abi + double ProtonEnergy = signal_proton.E() / 1000.0; // GeV Abi; + double ProtonKE = ProtonEnergy - ProtonMass_GeV; + + // ECal energy reconstruction + double ECal = ((MuonEnergy) + (ProtonKE) + BE); // GeV + // QE Energy Reconstruction + + double EQENum = 2 * (NeutronMass_GeV - BE) * MuonEnergy - + (BE * BE - 2 * NeutronMass_GeV * BE + + MuonMass_GeV * MuonMass_GeV + DeltaM2); + double EQEDen = + 2 * (NeutronMass_GeV - BE - MuonEnergy + vpmu.Mag() * vpmu.CosTheta()); + double EQE = EQENum / EQEDen; + + // 3D KI variables + + // For the calculation of the masses + // https://journals.aps.org/prc/pdf/10.1103/PhysRevC.95.065501 + + double MA = (22 * NeutronMass_GeV) + (18 * ProtonMass_GeV) - 0.34381; // GeV + + // For the calculation of the excitation energies + // https://doi.org/10.1140/epjc/s10052-019-6750-3 + + double MAPrime = + MA - NeutronMass_GeV + 0.0309; // GeV, constant obtained from table 7 + + // For the calculation of p_n, back to the Minerva PRL + // https://journals.aps.org/prl/pdf/10.1103/PhysRevLett.121.022504 + + double R = + MA + vpmuL.Mag() + vpL.Mag() - MuonEnergy - ProtonEnergy; // Equation 8 + + // Equation 7 + + double PL = + 0.5 * R - (MAPrime * MAPrime + vSumT.Mag() * vSumT.Mag()) / (2 * R); + + double DeltaPn = + TMath::Sqrt((vSumT.Mag() * vSumT.Mag()) + (PL * PL)) / 1000.0; + // https://journals.aps.org/prd/pdf/10.1103/PhysRevD.101.092001 + + TVector3 UnitZ(0, 0, 1); + double DeltaPtx = ( UnitZ.Cross(vpmuT) ).Dot(vSumT) / vpmuT.Mag() / 1000.0; // GeV/c + double DeltaPty = -(vpmuT).Dot(vSumT) / vpmuT.Mag() / 1000.0; // GeV/c + + //----------------------------------------// + + if (fDist == kDeltaPT) { + fXVar = DeltaPT; + } + else if (fDist == kDeltaAlphaT) { + fXVar = DeltaAlphaT; + } + else if (fDist == kDeltaPhiT) { + fXVar = DeltaPhiT; + } + else if (fDist == kMuonCosTheta) { + fXVar = MuonCosTheta; + } + else if (fDist == kProtonCosTheta) { + fXVar = ProtonCosTheta; + } + else if (fDist == kMuonMomentum) { + fXVar = MuonMomentum; + } + else if (fDist == kProtonMomentum) { + fXVar = ProtonMomentum; + } + else if (fDist == kDeltaPn) { + fXVar = DeltaPn; + } + else if (fDist == kDeltaPtx) { + fXVar = DeltaPtx; + } + else if (fDist == kDeltaPty) { + fXVar = DeltaPty; + } + else if (fDist == kECal) { + fXVar = ECal; + } + else if (fDist == kEQE) { + fXVar = EQE; + } + else { + NUIS_ABORT( + "MicroBooNE_CC1Mu1p_XSec_2D_nu: Didn't get a valid distribution"); + } + + if (fSlice == kDeltaPT) { + fYVar = DeltaPT; + } + else if (fSlice == kDeltaAlphaT) { + fYVar = DeltaAlphaT; + } + else if (fSlice == kDeltaPhiT) { + fYVar = DeltaPhiT; + } + else if (fSlice == kMuonCosTheta) { + fYVar = MuonCosTheta; + } + else if (fSlice == kProtonCosTheta) { + fYVar = ProtonCosTheta; + } + else if (fSlice == kMuonMomentum) { + fYVar = MuonMomentum; + } + else if (fSlice == kProtonMomentum) { + fYVar = ProtonMomentum; + } + else if (fSlice == kDeltaPn) { + fYVar = DeltaPn; + } + else if (fSlice == kDeltaPtx) { + fYVar = DeltaPtx; + } + else if (fSlice == kDeltaPty) { + fYVar = DeltaPty; + } + else if (fSlice == kECal) { + fYVar = ECal; + } + else if (fSlice == kEQE) { + fYVar = EQE; + } + else { + NUIS_ABORT( + "MicroBooNE_CC1Mu1p_XSec_2D_nu: Didn't get a valid slice"); + } + +} + +void MicroBooNE_CC1Mu1p_XSec_2D_nu::FillHistograms() { + Measurement1D::FillHistograms(); + if (Signal) { + fMCHist_Fine2D->Fill(fXVar, fYVar, Weight); + FillMCSlice(fXVar, fYVar, Weight); + } +} + +void MicroBooNE_CC1Mu1p_XSec_2D_nu::FillMCSlice(double x, double y, + double w) { + // Fill corresponding MC slice histogram + int slice_id = 0; + for (const auto& slice : fBinScheme.GetSlices()) { + std::pair edges = slice.GetSliceEdges(); + if (y >= edges.first && y < edges.second) { + fMCHist_Slices[slice_id]->Fill(x, w); + break; + } + slice_id++; + } +} + +void MicroBooNE_CC1Mu1p_XSec_2D_nu::ConvertEventRates() { + + // Do standard conversion + Measurement1D::ConvertEventRates(); + + // Apply MC truth -> reco smearing + std::vector slices_true; + for (size_t i = 0; i < fMCHist_Slices.size(); i++) { + TH1D *h = (TH1D *)fMCHist_Slices[i]->Clone(TString(fMCHist_Slices[i]->GetName()) + "_true"); + slices_true.push_back(h); + } + + for (int ireco = 1; ireco < fMCHist->GetNbinsX() + 1; ireco++) { + double total = 0.0; + for (int itrue = 1; itrue < fMCHist->GetNbinsX() + 1; itrue++) { + Slice slice = fBinScheme.GetSliceFromGlobal(itrue); + int slice_id = slice.GetSliceId(); + int itrue_local = itrue - slice.GetGlobalBins()[0] + 1; + TH1D *h = slices_true[slice_id]; + total += h->GetBinContent(itrue_local) * + h->GetBinWidth(itrue_local) * + fSmearingMatrix->operator()(ireco - 1, itrue - 1); + } + Slice slice = fBinScheme.GetSliceFromGlobal(ireco); + int slice_id = slice.GetSliceId(); + int ireco_local = ireco - slice.GetGlobalBins()[0] + 1; + TH1D *h = fMCHist_Slices[slice_id]; + h->SetBinContent(ireco_local, total / h->GetBinWidth(ireco_local)); + } + + for (size_t i = 0; i < slices_true.size(); i++) { + delete slices_true[i]; + } + + // Scale MC slices also by their width in Y + std::vector slices = fBinScheme.GetSlices(); + for (size_t i = 0; i < fMCHist_Slices.size(); i++) { + std::pair edges = slices[i].GetSliceEdges(); + float w = edges.second - edges.first; + fMCHist_Slices[i]->Scale(1.0 / w); + } + + // Convert into 1D list + fMCHist->Reset(); + int bincount = 0; + for (size_t i = 0; i < fDataHist_Slices.size(); i++) { + for (int j = 0; j < fDataHist_Slices[i]->GetNbinsX(); j++) { + fMCHist->SetBinContent(bincount + 1, fMCHist_Slices[i]->GetBinContent(j + 1)); + fMCHist->SetBinError(bincount + 1, fMCHist_Slices[i]->GetBinError(j + 1)); + bincount++; + } + } + +} + +void MicroBooNE_CC1Mu1p_XSec_2D_nu::SetHistograms() { + + // Open input file + std::string inputFileName = FitPar::GetDataBase() + + "/MicroBooNE/CC1Mu1p/CC1Mu1p_DataRelease.root"; + TFile *inputFile = TFile::Open(inputFileName.c_str()); + assert(inputFile && inputFile->IsOpen()); + + // Read bin configuration from file + LoadBinScheme(); + + std::string sample_name = fSettings.GetS("name"); + + // Get data histogram + fDataHist = (TH1D*)inputFile->Get(("TotalUnc_" + fSuffix).c_str()); + fDataHist->SetNameTitle(Form("%s_data", sample_name.c_str()), + Form("%s_data%s", sample_name.c_str(), fSettings.GetFullTitles().c_str())); + fDataHist->Scale(1e-38); + int nbins = fDataHist->GetNbinsX(); + + // Get covariance and smearing matrices + fFullCovar = new TMatrixDSym(fDataHist->GetNbinsX()); + fSmearingMatrix = new TMatrixD(fDataHist->GetNbinsX(), fDataHist->GetNbinsX()); + // Convert from TH2D to TMatrixD + TH2D* temp_cov = (TH2D*)inputFile->Get(("Cov_" + fSuffix).c_str()); + TH2D* temp_smr = (TH2D*)inputFile->Get(("Ac_" + fSuffix).c_str()); + for (int i = 0; i < nbins; ++i) { + for (int j = 0; j < nbins; ++j) { + (*fFullCovar)(i, j) = temp_cov->GetBinContent(i + 1, j + 1); + (*fSmearingMatrix)(i, j) = temp_smr->GetBinContent(i + 1, j + 1); + } + } + + // Set fine-grained MC 2D hist + double slicesMin = fBinScheme.GetSlices().front().GetSliceEdges().first; + double slicesMax = fBinScheme.GetSlices().back().GetSliceEdges().second; + double binsMin = DBL_MAX, binsMax = -DBL_MAX; + for (auto& s : fBinScheme.GetSlices()) { + binsMin = std::min(binsMin, s.GetBinsForTH1()[0]); + binsMax = std::max(binsMax, s.GetBinsForTH1()[s.GetNumberBins()]); + } + fMCHist_Fine2D = new TH2D(Form("%s_MC_FINE_2D", sample_name.c_str()), + Form("%s_MC_FINE_2D; %s; %s; %s", + sample_name.c_str(), + fSettings.GetXTitle().c_str(), fSliceTitle.c_str(), fSettings.GetYTitle().c_str()), + 100, binsMin, binsMax, 100, slicesMin, slicesMax); + SetAutoProcessTH1(fMCHist_Fine2D); + + // Divide data in slices + int slice_id = 0; + for (const auto& slice : fBinScheme.GetSlices()) { + + // Add data histogram to slice + TString name = Form("%s_data_Slice%lu", sample_name.c_str(), slice_id); + TString title = Form("%s_data_Slice%lu; %s; %s", + sample_name.c_str(), slice_id, + fSettings.GetXTitle().c_str(), fSettings.GetYTitle().c_str()); + TH1D* temp_data = new TH1D(name, title, slice.GetNumberBins(), slice.GetBinsForTH1()); + + std::vector temp_global_bins = slice.GetGlobalBins(); + double slice_width = slice.GetSliceWidth(); + for (int i = 1; i <= temp_global_bins.size(); ++i) { + // Get bin properties + int global_bin = temp_global_bins[i-1]; + double content = fDataHist->GetBinContent(global_bin); + double error = fDataHist->GetBinError(global_bin); + double width = temp_data->GetBinWidth(i); // width in slice hist only! + // Compute new value + double scale_factor = 1.0 / width / slice_width; + double new_content = content * scale_factor; + double new_error = error * scale_factor; + // Fill slice histogram + temp_data->SetBinContent(i, new_content); + temp_data->SetBinError(i, new_error); + // Update global data histogram + fDataHist->SetBinContent(global_bin, new_content); + fDataHist->SetBinError(global_bin, new_error); + } + + temp_data->Sumw2(); + fDataHist_Slices.push_back(temp_data); + + // Add MC histogram to slice + fMCHist_Slices.push_back((TH1D *)temp_data->Clone()); + name = Form("%s_MC_Slice%lu", sample_name.c_str(), slice_id); + title = Form("%s_MC_Slice%lu; %s; %s", + sample_name.c_str(), slice_id, + fSettings.GetXTitle().c_str(), fSettings.GetYTitle().c_str()); + fMCHist_Slices[slice_id]->SetNameTitle(name, title); + fMCHist_Slices[slice_id]->Reset(); + + SetAutoProcessTH1(fDataHist_Slices[slice_id], kCMD_Write); + SetAutoProcessTH1(fMCHist_Slices[slice_id]); + + slice_id++; + + } + + // Update covariance matrix + for (int i = 0; i < fBinScheme.GetNumberBins(); ++i) { + Slice iSlice = fBinScheme.GetSliceFromGlobal(i+1); + int iLocal = i + 1 - iSlice.GetGlobalBins()[0] + 1; + double iBinWidth = fDataHist_Slices[iSlice.GetSliceId()]->GetBinWidth(iLocal); + double iSliceWidth = iSlice.GetSliceWidth(); + for (int j = 0; j < fBinScheme.GetNumberBins(); ++j) { + Slice jSlice = fBinScheme.GetSliceFromGlobal(j+1); + int jLocal = j + 1 - jSlice.GetGlobalBins()[0] + 1; + double jBinWidth = fDataHist_Slices[jSlice.GetSliceId()]->GetBinWidth(jLocal); + double jSliceWidth = jSlice.GetSliceWidth(); + double cov_scale_factor = 1.0 / iBinWidth / jBinWidth / iSliceWidth / jSliceWidth; + (*fFullCovar)(i, j) = (*fFullCovar)(i, j) * cov_scale_factor; + } + } + + fDecomp = StatUtils::GetDecomp(fFullCovar); + +} + +void MicroBooNE_CC1Mu1p_XSec_2D_nu::LoadBinScheme() { + + // Open input file + std::string binFileName = FitPar::GetDataBase() + + "/MicroBooNE/CC1Mu1p/CC1Mu1p_BinScheme.txt"; + std::ifstream binFile(binFileName); + + std::string line; + bool inTargetBlock = false; + + while (std::getline(binFile, line)) { + + // Trim leading/trailing spaces + line.erase(0, line.find_first_not_of(" \t")); + line.erase(line.find_last_not_of(" \t\r\n") + 1); + + // Detect start of block + if (line.find(fSuffix+"Plot") != std::string::npos) { + inTargetBlock = true; + continue; + } + + // Stop when next block begins + if (inTargetBlock && line.rfind("Serial", 0) == 0 && line != fSuffix+"Plot") { + break; + } + + // Skip header or empty lines + if (!inTargetBlock || line.empty() || line[0] == 'B') continue; + + // Split line by '&' + std::stringstream ss(line); + std::string token; + std::vector parts; + while (std::getline(ss, token, '&')) { + // Trim whitespace around each part + token.erase(0, token.find_first_not_of(" \t")); + token.erase(token.find_last_not_of(" \t\r\n") + 1); + parts.push_back(token); + } + + // Add bin to bin scheme + fBinScheme.AddBin(parts); + + } + + binFile.close(); + +} \ No newline at end of file diff --git a/src/MicroBooNE/MicroBooNE_CC1Mu1p_XSec_2D_nu.h b/src/MicroBooNE/MicroBooNE_CC1Mu1p_XSec_2D_nu.h new file mode 100644 index 000000000..bcd369fb8 --- /dev/null +++ b/src/MicroBooNE/MicroBooNE_CC1Mu1p_XSec_2D_nu.h @@ -0,0 +1,140 @@ +// Copyright 2016 L. Pickering, P Stowell, R. Terri, C. Wilkinson, C. Wret + +/******************************************************************************* +* This file is part of NUISANCE. +* +* NUISANCE is free software: you can redistribute it and/or modify +* it under the terms of the GNU General Public License as published by +* the Free Software Foundation, either version 3 of the License, or +* (at your option) any later version. +* +* NUISANCE is distributed in the hope that it will be useful, +* but WITHOUT ANY WARRANTY; without even the implied warranty of +* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +* GNU General Public License for more details. +* +* You should have received a copy of the GNU General Public License +* along with NUISANCE. If not, see . +*******************************************************************************/ +#ifndef MICROBOONE_CC1MU1P_2D_NU_H_SEEN +#define MICROBOONE_CC1MU1P_2D_NU_H_SEEN + +#include +#include "Measurement1D.h" + +#include +#include +#include + +class TH2D; + +class Slice { +public: + // Constructor + Slice(std::string name, int id); + // Virtual destructor + ~Slice(); + // Add single bin to slice + void AddBin(int global, double low, double high); + // Get slice name + inline std::string GetSliceName() const { return slice_name; }; + // Get slice index + inline int GetSliceId() const { return slice_id; }; + // Get slice edges + inline std::pair GetSliceEdges() const { return slice_edges; }; + inline double GetSliceWidth() const { return slice_edges.second - slice_edges.first; }; + // Get global bin numbers for slice + inline int GetNumberBins() const { return global_bins.size(); }; + inline std::vector GetGlobalBins() const { return global_bins; }; + // Get bin edges for slice (different formats) + inline std::vector> GetBins() const { return bin_edges; }; + double* GetBinsForTH1() const; + +private: + std::string slice_name; + int slice_id; + std::pair slice_edges; + std::vector> bin_edges; + std::vector global_bins; + +}; + +class BinScheme { +public: + // Constructor + BinScheme(); + // Virtual destructor + ~BinScheme(); + // Add single bin to scheme + void AddBin(std::vector bin_config); + void AddBin(std::string slice_name, int global, double low, double high); + // Get all slices + inline std::vector GetSlices() { return slice_vec; }; + // Get total number of bins + int GetNumberBins(); + // Get slice containing given global bin number + Slice GetSliceFromGlobal(int global); + +private: + int current_slice = 0; + std::vector slice_vec; +}; + +class MicroBooNE_CC1Mu1p_XSec_2D_nu : public Measurement1D { +public: + /// Basic Constructor + MicroBooNE_CC1Mu1p_XSec_2D_nu(nuiskey samplekey); + + /// Virtual Destructor + ~MicroBooNE_CC1Mu1p_XSec_2D_nu() {}; + + /// Apply signal definition + bool isSignal(FitEvent* nvect); + + /// Read histograms in a special way because format is different + void SetHistograms(); + + /// Fill kinematic distributions + void FillEventVariables(FitEvent* customEvent); + + // Read data histograms from file + void FillHistograms(); + + /// Additional smearing matrix multiplication by Ac + void ConvertEventRates(); + + /// Read bin configuration from file + void LoadBinScheme(); + +private: + enum Distribution { kDeltaPT=0, + kDeltaAlphaT=1, + kDeltaPhiT=2, + kMuonCosTheta=3, + kProtonCosTheta=4, + kMuonMomentum=5, + kProtonMomentum=6, + kDeltaPn=7, + kDeltaPtx=8, + kDeltaPty=9, + kECal=10, + kEQE=11 }; + + std::vector fMCHist_Slices; + std::vector fDataHist_Slices; + TH2D *fMCHist_Fine2D; + TMatrixD* fSmearingMatrix; + + // Input distribution information + std::string fSuffix; + Distribution fSlice; + Distribution fDist; + std::string fSliceTitle; + + BinScheme fBinScheme; + + void FillMCSlice(double x, double y, double w); + +}; + +#endif