root-project / root

The official repository for ROOT: analyzing, storing and visualizing big data, scientifically
https://root.cern
Other
2.63k stars 1.25k forks source link

[math] Minuit2 errors out with `Initial matrix not pos.def.` depending on which histograms were fitted before #13852

Closed guitargeek closed 10 months ago

guitargeek commented 11 months ago

This is the underlying issue related to a CMSSW issue:

To reproduce, download this ROOT file with a single TH2F: histo.zip

And then fit a subset of the 1D slices:

std::unique_ptr<TFile> file{TFile::Open("histo.root", "READ")};
auto hist = file->Get<TH2>("dxyres_vs_eta");
hist->FitSlicesY(nullptr, 10, 21, 0, "QNR");

You will see a Minuit2 error:

Error in <Minuit2>: VariableMetricBuilder Initial matrix not pos.def.

Note that this doesn't happen with the old Minuit, which you can select like this at the beginning of the script:

ROOT::Math::MinimizerOptions::SetDefaultMinimizer("Minuit", "Migrad");

There is something weird going on, because this happens in the fit for bin 21, but it is necessary to fit all the slices 10 to 20 before. Otherwise there is no error. In other words, fitting the TH1 from the 21th slice in isolation doesn't result in the error.

lmoneta commented 10 months ago

If you try fitting the obtained slice (slice 21), you will see that in both Minuit and Minuit2 the fit did not work. Just run this simple code:

   auto f = TFile::Open("histo.root", "READ");
   auto hist = f->Get<TH2>("dxyres_vs_eta");

   auto h20 = hist->ProjectionY("h20",20,20);
   auto h21 = hist->ProjectionY("h21",21,21);
   auto c1 = new TCanvas();
   c1->Divide(1,2);
   c1->cd(1);
   h20->Fit("gaus");
   c1->cd(2);
   // second fit fails                                                                                                                                              
   h21->Fit("gaus");

If you run only the second fit, it works because some default steps sizes are used at the beginning. You will get better slice fits if using option L when fitting the slices:

hist->FitSlicesY(nullptr, 10, 21, 0, "LR");

and defining a restricted range for the fitted functions to avoid fitting the outlier events

Close the issue since it is not a bug.