Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
31 changes: 22 additions & 9 deletions CIValidations/SplinePlotting/SplinePlotting.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -33,16 +33,28 @@ TSpline3* MakeSpline(double x[], double y[], int N,
void MakeSplinePlot(double x[], double y[], int N, TCanvas* c) {
constexpr int nInterpolations = 5;
constexpr SplineInterpolation interpolationTypes[nInterpolations] = {kTSpline3, kMonotonic, kLinear, kAkima, kKochanekBartels};
constexpr Color_t colors[nInterpolations] = {kBlue, kRed, kGreen, kViolet, kCyan};
constexpr Color_t colors[nInterpolations] = {kBlue, kRed, kGreen, kViolet, kAzure};
constexpr Style_t lineStyles[nInterpolations] = {kSolid, kDotted, kDashed, kDashDotted, kDotted};

std::vector<TSpline3*> splines;

// Create and draw TGraph to show original points (knots)
TGraph* knotsGraph = new TGraph(N, x, y);
auto knotsGraph = std::make_unique<TGraph>(N, x, y);
knotsGraph->SetMarkerStyle(kFullCircle);
knotsGraph->SetMarkerSize(1.2);
knotsGraph->SetMarkerColor(kBlack);
// Set titles for the axes
knotsGraph->GetXaxis()->SetTitle("Value");
knotsGraph->GetYaxis()->SetTitle("Weight");

// Find the minimum and maximum y-values in the data
double y_min = *std::min_element(y, y + N);
double y_max = *std::max_element(y, y + N);

// Add 10% padding to the y-axis range
double padding = 0.20 * (y_max - y_min);
knotsGraph->GetYaxis()->SetRangeUser(y_min - padding, y_max + padding);

knotsGraph->Draw("AP"); // "A" creates new axis, "P" draws points

// Draw each spline
Expand All @@ -57,8 +69,8 @@ void MakeSplinePlot(double x[], double y[], int N, TCanvas* c) {
auto legend = std::make_unique<TLegend>(0.15, 0.75, 0.45, 0.88);
legend->SetFillStyle(0);
legend->SetBorderSize(0);
legend->SetTextSize(0.03);
legend->AddEntry(knotsGraph, "Knots", "p");
legend->SetTextSize(0.04);
legend->AddEntry(knotsGraph.get(), "Knots", "p");

for (int i = 0; i < nInterpolations; ++i) {
legend->AddEntry(splines[i], SplineInterpolation_ToString(interpolationTypes[i]).c_str(), "l");
Expand All @@ -70,7 +82,6 @@ void MakeSplinePlot(double x[], double y[], int N, TCanvas* c) {
for (int i = 0; i < nInterpolations; ++i) {
delete splines[i];
}
delete knotsGraph;
}


Expand Down Expand Up @@ -124,10 +135,13 @@ void MakeTF1Plot(double x[], double y[], int N, TCanvas* c) {
std::vector<TF1*> splines;

// Create and draw TGraph to show original points (knots)
TGraph* knotsGraph = new TGraph(N, x, y);
auto knotsGraph = std::make_unique<TGraph>(N, x, y);
knotsGraph->SetMarkerStyle(kFullCircle);
knotsGraph->SetMarkerSize(1.2);
knotsGraph->SetMarkerColor(kBlack);
// Set titles for the axes
knotsGraph->GetXaxis()->SetTitle("Value");
knotsGraph->GetYaxis()->SetTitle("Weight");
knotsGraph->Draw("AP"); // "A" creates new axis, "P" draws points

// Draw each spline
Expand All @@ -140,10 +154,10 @@ void MakeTF1Plot(double x[], double y[], int N, TCanvas* c) {

// Add legend
auto legend = std::make_unique<TLegend>(0.15, 0.75, 0.45, 0.88);
legend->SetTextSize(0.03);
legend->SetTextSize(0.04);
legend->SetFillStyle(0);
legend->SetBorderSize(0);
legend->AddEntry(knotsGraph, "Knots", "p");
legend->AddEntry(knotsGraph.get(), "Knots", "p");

for (int i = 0; i < nInterpolations; ++i) {
legend->AddEntry(splines[i],
Expand All @@ -156,7 +170,6 @@ void MakeTF1Plot(double x[], double y[], int N, TCanvas* c) {
for (int i = 0; i < nInterpolations; ++i) {
delete splines[i];
}
delete knotsGraph;
}

void TF1_Plot(TCanvas* c) {
Expand Down
39 changes: 39 additions & 0 deletions CIValidations/UnitTests/histogram_tests.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -156,3 +156,42 @@ TEST_CASE("TH2Poly integral", "[HistogramUtils]") {

delete Poly;
}

// KS: I have trust issues with respect to TH2Poly
TEST_CASE("TH2Poly Bin Geometry", "[HistogramUtils]") {
TH2Poly* Poly = GetPolyHist();

// Check bin 10 (ROOT bins are 1-based, TH2Poly internal bins are 0-based)
const int bin = 10;

// ROOT quirk: GetBinContent(10) corresponds to GetBins(9)
TH2PolyBin* Bin = static_cast<TH2PolyBin*>(Poly->GetBins()->At(bin - 1));

// Content must match
REQUIRE_THAT(Poly->GetBinContent(bin),
Catch::Matchers::WithinAbs(Bin->GetContent(), 1e-12));

// Extract geometry
const double xmin = Bin->GetXMin();
const double xmax = Bin->GetXMax();
const double ymin = Bin->GetYMin();
const double ymax = Bin->GetYMax();

// Get the original binning definition
std::vector<double> BinningX;
std::vector<double> BinningY;
GetBinning(BinningX, BinningY);

// Bin numbering in GetPolyHist():
// loops are iy (outer) then ix (inner), so:
// globalBin = iy*(Nx-1) + ix + 1
const int ix = (bin - 1) % (BinningX.size() - 1);
const int iy = (bin - 1) / (BinningX.size() - 1);

REQUIRE_THAT(xmin, Catch::Matchers::WithinAbs(BinningX[ix], 1e-12));
REQUIRE_THAT(xmax, Catch::Matchers::WithinAbs(BinningX[ix + 1], 1e-12));
REQUIRE_THAT(ymin, Catch::Matchers::WithinAbs(BinningY[iy], 1e-12));
REQUIRE_THAT(ymax, Catch::Matchers::WithinAbs(BinningY[iy + 1], 1e-12));

delete Poly;
}
73 changes: 52 additions & 21 deletions README.md
Original file line number Diff line number Diff line change
Expand Up @@ -8,9 +8,8 @@ This tutorial will help you get used to different parts of MaCh3.

MaCh3 is predominantly C++ software, although some functionality is available through python as well. See the table of contents for more.

[![Code - Documented](https://img.shields.io/badge/Code-Documented-2ea44f)](https://github.com/mach3-software/MaCh3/wiki)
[![Code - Documented](https://img.shields.io/badge/Code-Documented-2ea44f)](https://mach3-software.github.io/MaCh3/index.html)
[![Container Image](https://img.shields.io/badge/Container-Image-brightgreen)](https://github.com/mach3-software/MaCh3Tutorial/pkgs/container/mach3tutorial)
[![Code - Doxygen](https://img.shields.io/badge/Code-Doxygen-2ea44f)](https://mach3-software.github.io/MaCh3/index.html)

## Table of contents
1. [How to Start?](#how-to-start)
Expand All @@ -31,7 +30,6 @@ MaCh3 is predominantly C++ software, although some functionality is available th
3. [Changing Oscillation Engine](#changing-oscillation-engine)
4. [Atmospheric Sample](#atmospheric-sample)
5. [Plotting Kinematic Distribution](#plotting-kinematic-distribution)
6. [More Advanced Development](#more-advanced-development)
6. [MCMC Diagnostic](#mcmc-diagnostic)
1. [Running Multiple Chains](#running-multiple-chains)
7. [Useful Settings](#useful-settings)
Expand All @@ -55,9 +53,9 @@ source bin/setup.MaCh3.sh
source bin/setup.MaCh3Tutorial.sh
```

If this does not work, check out the [Requirements](https://github.com/mach3-software/MaCh3/tree/Charlotte-Knight/tidy/MaCh3-on-Mac?tab=readme-ov-file#system-requirements) for MaCh3 in case you are missing something.
If this does not work, check out the [Requirements](https://github.com/mach3-software/MaCh3?tab=readme-ov-file#system-requirements) for MaCh3 in case you are missing something.

One of the suggestions in [Requirements](https://github.com/mach3-software/MaCh3/tree/Charlotte-Knight/tidy/MaCh3-on-Mac?tab=readme-ov-file#system-requirements) is to use a conda/micromamba environment. You can install everything you need and build the tutorial using these commands:
One of the suggestions in [Requirements](https://github.com/mach3-software/MaCh3?tab=readme-ov-file#system-requirements) is to use a conda/micromamba environment. You can install everything you need and build the tutorial using these commands:

<details>
<summary>Click to see commands</summary>
Expand Down Expand Up @@ -110,7 +108,7 @@ Alternatively you can use containers by
```bash
docker pull ghcr.io/mach3-software/mach3tutorial:alma9latest
```
To read more about how to use containers, check our wiki [here](https://github.com/mach3-software/MaCh3/wiki/12.-Containers).
To read more about how to use containers, check our doxygen [here](https://mach3-software.github.io/MaCh3/Containers.html).

## How to run MCMC
To run MCMC simply
Expand Down Expand Up @@ -156,7 +154,7 @@ where `Test.root` is the output of running `MCMCTutorial` as described [here](#h

This is the marginalised posterior of a single parameter. This is the main output of the MCMC.

There are many options in `ProcessMCMC` that allow you to make many more analysis plots from the MCMC output; we recommend that you see [here](https://github.com/mach3-software/MaCh3/wiki/09.-Bayesian-Analysis,-Plotting-and-MCMC-Processor) to get better idea what each plot means. In particular, we recommend comparing 2D posteriors with the correlation matrix, and playing with triangle plots.
There are many options in `ProcessMCMC` that allow you to make many more analysis plots from the MCMC output; we recommend that you see [here](https://mach3-software.github.io/MaCh3/BayesianAnalysis.html) to get better idea what each plot means. In particular, we recommend comparing 2D posteriors with the correlation matrix, and playing with triangle plots.

### Plotting MCMC Posteriors
Once you have the processed MCMC output from `ProcessMCMC`, which will be called something like `<inputName>_Process.root`, you can make fancier analysis plots from it using the `GetPostfitParamPlots` app like:
Expand All @@ -177,7 +175,7 @@ The output should look like plot below. This conveys same information as the ind

<img width="350" alt="Ridge" src="https://github.com/user-attachments/assets/617f5929-b389-495e-ab7b-2ecd6c2d991e">

**Violin Plot** - This also allow to see nicely non-Gaussian parameters but also is usefull in comparing two chains.
**Violin Plot** - This also allow to see nicely non-Gaussian parameters but also is useful in comparing two chains.
`ProcessMCMC` must be run with option "PlotCorr" to be able to produce violin plot.

<img width="350" alt="Violin example" src="https://github.com/user-attachments/assets/4788ab29-f24a-4b09-8b0f-c9b36d069cfe">
Expand Down Expand Up @@ -288,11 +286,11 @@ First let's get a better understanding `TutorialConfigs/CovObjs/SystematicModel.
MCMC: 0.2
Type: Norm
```
If you want to read more about the specifics of implementating systematics in the configs, please see [here](https://github.com/mach3-software/MaCh3/wiki/02.-Implementation-of-Systematic)
If you want to read more about the specifics of implementing systematics in the configs, please see [here](https://mach3-software.github.io/MaCh3/group__SamplesAndParameters.html)

As a first step let's modify `Error: 0.11` to `Error: 2.0`. This significantly changes the prior uncertainty on the `Norm_Param_0` parameter. Such a change should be very noticeable in the MCMC posterior.

Let's also modify the output file name so we don't overwrite our previous outputs. This is governed by the manager class (read more [here](https://github.com/mach3-software/MaCh3/wiki/01.-Manager-and-config-handling)).
Let's also modify the output file name so we don't overwrite our previous outputs. This is governed by the manager class.
You can modify this in the configs by for example changing `OutputFile: "Test.root"` in `TutorialConfigs/FitterConfig.yaml` to `OutputFile: "Test_Modified.root"` (or you could use a parameter override by adding `General:OutputFile:Test_Modified.root` as a command line argument).

With the modifications, let's run MCMC again:
Expand Down Expand Up @@ -340,15 +338,43 @@ SelectionCuts:

You can also easily change what variable the sample is binned in to get slightly different fit results. For example, we could change `RecoNeutrinoEnergy` to `TrueNeutrinoEnergy` in `TutorialConfigs/Samples/SampleHandler_Tutorial.yaml`:
```yaml
Binning:
XVarStr : "TrueNeutrinoEnergy"
XVarBins: [0., 0.5, 1., 1.25, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.25, 3.5, 3.75, 4., 5., 6., 10.]
Binning:
VarStr : [ "RecoNeutrinoEnergy" ]
VarBins: [[0., 0.5, 1., 1.25, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.25, 3.5, 3.75, 4., 5., 6., 10.]]
Uniform: true
```
You can run the MCMC again using this new configuration (after changing the output file name again!), and then compare all 3 chains using:
```bash
./bin/ProcessMCMC bin/TutorialDiagConfig.yaml Test.root Default_Chain Test_Modified.root Modified_Chain Test_Modified_Sample.root ModifiedSameple_Chain
```

<details>
<summary><strong>(Detailed) Uniform vs NonUniform binning</strong></summary>
Uniform binning, discussed above creating kind of uniform grid, can be easily scaled like this
```yaml
Binning:
VarStr : ["RecoNeutrinoEnergy", "TrueQ2"]
VarBins: [ [0., 0.5, 1., 1.25, 1.5, 1.75, 2., 2.25, 2.5, 2.75, 3., 3.25, 3.5, 3.75, 4., 5., 6., 10.],
[0., 0.5, 1., 1.25, 1.5, 1.75, 2., 5] ]
Uniform: true
```

However there is option to use non-unform binning which characterises by having bin sizes vary,
but all bins are axis-aligned hyper-rectangles.
We recommend to read more [here](https://mach3-software.github.io/MaCh3/classBinningHandler.html).

Example of such binning is given below:
```yaml
Binning:
VarStr : ["RecoNeutrinoEnergy", "TrueQ2"]
Bins: [[[0.0, 0.1], [0.0, 0.1]], [[0.0, 0.1], [0.1, 0.2]], [[0.0, 0.1], [0.2, 0.3]],
[[0.0, 0.1], [0.3, 0.4]], [[0.0, 0.1], [0.4, 0.5]], [[0.0, 0.1], [0.5, 0.6]],
[[0.0, 0.1], [0.6, 0.7]], [[0.0, 0.1], [0.7, 0.8]], [[0.0, 0.1], [0.8, 0.9]],
[[0.0, 0.1], [0.9, 1.0]] ]
Uniform: false
```
</details>

### Adding a New Sample
Up to this point we have only modified an existing sample, but how would we add a new one? We can start by first making a copy of the sample config `TutorialConfigs/Samples/SampleHandler_Tutorial.yaml` and calling it `TutorialConfigs/Samples/SampleHandler_User.yaml`.
For the moment feel free to change the sample name, binning, whatever you want. Go wild! Just be sure to keep `TutorialConfigs/CovObjs` the same, or things will break.
Expand Down Expand Up @@ -382,10 +408,13 @@ This approach may give a performance boost, but is especially nice if you are sh
### Changing Oscillation Engine
MaCh3 has access to many oscillation engines via the `NuOscillator` framework. You can check the features of this using:
```bash
bin/mach3-config --features
mach3-config --features
MULTITHREAD MR2T2 PSO Minuit2 Prob3ppLinear NuFast
```
Thus, you can easily access information about MaCh3 features, most importantly the fitter engines (`MR2T2`, `PSO`, `Minuit2`) and neutrino oscillation engines (`Prob2ppLinear`, `NuFast`).
Thus, you can easily access information about MaCh3 features, most importantly the fitter engines (`MR2T2`, `PSO`, `Minuit2`) and neutrino oscillation engines (`Prob3ppLinear`, `NuFast`).

> [!NOTE]
> To see all available engines that MaCh3 can access can be found [here](https://github.com/mach3-software/MaCh3/tree/develop?tab=readme-ov-file#oscillator)

By default, the oscillation engine we use is `NuFast`. However, you can change to another engine (`Prob3++` here) by modifying the sample config `TutorialConfigs/Samples/SampleHandler_Tutorial.yaml`:
```yaml
Expand Down Expand Up @@ -420,17 +449,19 @@ An example of what you might expect for output can be seen here:

<img width="350" alt="Kinematic example" src="https://github.com/user-attachments/assets/534bcb17-f26c-4fc2-a77a-5d253b0ed241">

### More Advanced Development
<details>
<summary><strong>(Detailed) More Advanced Development</strong></summary>
Not everything can be done by modifying config in sample implementation. The actual implementation of samples is in `Samples/SampleHandlerTutorial`, which inherits from `SampleHandlerFDBase`.
The latter class deals with actual event reweighting and general heavy lifting. `SampleHandlerTutorial` deals with more specific information, like MC variable loading. This is because each experiment has slightly different MC format and different information available, and so it must be implemented in a custom manner.
</details>

## MCMC Diagnostic
A crucial part of MCMC is diagnosing whether a chain has converged or not. You can produce chain diagnostics by running:
```bash
./bin/DiagMCMC Test.root bin/TutorialDiagConfig.yaml
```
This will produce a plethora of diagnostic information. One of most important of these are the autocorrelations for each of the parameters. Autocorrelation indicates how correlated MCMC steps are when they are n-steps apart in the chain.
Generally speaking, we want the autocorrelation to drop to 0 fast and stay around there (like in the plot below). You can read about other diagnostics [here](https://github.com/mach3-software/MaCh3/wiki/11.-Step-size-tuning), but it is sufficient for now to focus on autocorrelation.
Generally speaking, we want the autocorrelation to drop to 0 fast and stay around there (like in the plot below). You can read about other diagnostics [here](https://mach3-software.github.io/MaCh3/MCMCconvergance.html), but it is sufficient for now to focus on autocorrelation.

<img width="350" alt="Posterior example" src="https://github.com/user-attachments/assets/e146698c-713c-4daf-90df-adeb3051539b">

Expand Down Expand Up @@ -461,7 +492,7 @@ You can make plots from the diagnostic output using:
```bash
./bin/PlotMCMCDiag Test_MCMC_Diag.root "mcmc_diagnostics"
```
If you add a second/third arguemnt, the macro can compare several files:
If you add a second/third argument, the macro can compare several files:
```bash
./bin/PlotMCMCDiag Test_MCMC_Diag.root "first file label" SecondFile_MCMC_Diag.root "second file label"
```
Expand Down Expand Up @@ -512,9 +543,9 @@ Read more [here](https://mach3-software.github.io/MaCh3/Structs_8h.html#a960da89

## Other Useful Plots

There are a number of apps included to make plots from the results of your fits, llh scans etc. You can find more details on them and how they work in the main MaCh3 wiki [here](https://github.com/mach3-software/MaCh3/wiki). There you will also find some instructions on how you can write your own plotting scripts.
There are a number of apps included to make plots from the results of your fits, llh scans etc. You can find more details on them and how they work in the main MaCh3 doxygen page [here](https://mach3-software.github.io/MaCh3/index.html). There you will also find some instructions on how you can write your own plotting scripts.

The plotting library is configured using yaml files. You can see some examples of such config files in the plotting directory, and a detailed explanation of them is given in [the wiki](https://github.com/mach3-software/MaCh3/wiki).
The plotting library is configured using yaml files. You can see some examples of such config files in the plotting directory, and a detailed explanation of them is given in [here](https://mach3-software.github.io/MaCh3/group__MaCh3Plotting.html).

Some examples on how to make some "standard" plots are given below.

Expand Down Expand Up @@ -559,7 +590,7 @@ If you have installed the python interface for MaCh3 as described
[here](https://github.com/mach3-software/MaCh3?tab=readme-ov-file#python)
then you can also use the provided python plotting module. The details on how
to write custom python scripts using this plotting module are detailed in
[the wiki](https://github.com/mach3-software/MaCh3/wiki/15.-Plotting#custom-plotting-scripts).
[here](https://mach3-software.github.io/MaCh3/group__MaCh3Plotting.html).
Here we will walk you through some example scripts.

For the examples here, we will use matplotlib and numpy. These can be installed using the provided [requirements.txt](requirements.txt) by doing
Expand Down
6 changes: 5 additions & 1 deletion cmake/Modules/CIValidations.cmake
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,11 @@ install(DIRECTORY
if(MaCh3Tutorial_UNITTESTS_ENABLED OR MaCh3Tutorial_Benchmark_ENABLED)
find_package(Catch2 QUIET)
if(NOT Catch2_FOUND)
CPMAddPackage("gh:catchorg/Catch2@3.5.2")
CPMAddPackage(
NAME Catch2
GITHUB_REPOSITORY catchorg/Catch2
VERSION 3.7.1
)
LIST(APPEND CMAKE_MODULE_PATH ${Catch2_SOURCE_DIR}/extras)
endif()

Expand Down