#include <NuMatrixFitterMinuit.h>
Public Member Functions | |
| NuMatrixFitterMinuit () | |
| NuMatrixFitterMinuit (const Double_t FluxPoT, const string ndData, const string fdData, const string helperFile, const string xsecFile) | |
| virtual | ~NuMatrixFitterMinuit () |
| virtual double | operator() (const std::vector< double > &x) const |
| virtual double | Up () const |
| virtual const Double_t | CalculateChi2 (const TH1D *fdPred, const TH1D *fdData) const |
| virtual const Double_t | CalculateLikelihood (const TH1D *fdPred, const TH1D *fdData) const |
| virtual void | CPTFit (const string inputFilename, const string fdDataFilename, const string outputFilename) |
| virtual NuMatrixOutput * | DoCPTFit (NuMatrixInput *mmInput, TH1D *numufdData, TH1D *numubarfdData) |
| virtual void | CPTAnalysis () |
| virtual void | SingleCPTPrediction (Double_t dm2, Double_t sn2, Double_t dm2bar, Double_t sn2bar) |
Private Member Functions | |
| virtual void | SetupPars () |
Private Attributes | |
| NuMatrixInput * | fInput |
| TFitterMinuit * | fFitter |
| TH1D * | fnumufdData |
| TH1D * | fnumubarfdData |
| Double_t | ffluxPoT |
| std::string | fndDataFilename |
| std::string | ffdDataFilename |
| std::string | fhelperFilename |
| std::string | fxsecFilename |
| std::string | foutputFilename |
|
|
Definition at line 47 of file NuMatrixFitterMinuit.cxx. 00048 {
00049 fndDataFilename = "";
00050 ffdDataFilename = "";
00051 fhelperFilename = "";
00052 fxsecFilename = "";
00053 //fmmNuMu = 0;
00054 //fmmNuMuBar = 0;
00055 fInput = 0;
00056
00057 fnumufdData = 0;
00058 fnumubarfdData = 0;
00059
00060 fFitter = new TFitterMinuit();
00061 fFitter->CreateMinimizer();
00062 fFitter->SetMinuitFCN(this);
00063 }
|
|
||||||||||||||||||||||||
|
Definition at line 66 of file NuMatrixFitterMinuit.cxx. References ffdDataFilename, fFitter, ffluxPoT, fhelperFilename, fInput, fndDataFilename, foutputFilename, and fxsecFilename. 00071 {
00072 ffluxPoT = FluxPoT;
00073 fndDataFilename = ndData;
00074 ffdDataFilename = fdData;
00075 fhelperFilename = helperFile;
00076 fxsecFilename = xsecFile;
00077 foutputFilename =
00078 "$SRT_PRIVATE_CONTEXT/NtupleUtils/results/MMOutputMinuit.root";
00079 //fmmNuMu = 0;
00080 //fmmNuMuBar = 0;
00081 fInput = 0;
00082
00083 fFitter = new TFitterMinuit();
00084 fFitter->CreateMinimizer();
00085 fFitter->SetMinuitFCN(this);
00086 }
|
|
|
Definition at line 89 of file NuMatrixFitterMinuit.cxx. 00090 {
00091 }
|
|
||||||||||||
|
Definition at line 116 of file NuMatrixFitterMinuit.cxx. 00118 {
00119 //Return Chi2.
00120
00121 Double_t chi2 = 0;
00122
00123 for (Int_t i=1; i<=fdPred->GetNbinsX(); ++i){
00124 //Bizarre limits because root histograms are silly
00125 Double_t mnu = fdPred->GetBinContent(i);
00126 Double_t dnu = fdData->GetBinContent(i);
00127 // cout << "Bin " << i
00128 // << ": MC " << mnu
00129 // << "; data " << dnu
00130 // << endl;
00131 if (mnu<0.0001){mnu = 0.0001;}
00132 chi2 += (dnu-mnu)*(dnu-mnu)/mnu;
00133 }
00134 return chi2;
00135 }
|
|
||||||||||||
|
Definition at line 138 of file NuMatrixFitterMinuit.cxx. Referenced by operator()(), and SingleCPTPrediction(). 00141 {
00142 //Aim to minimise -2lnL. That's what I'm returning.
00143
00144 Double_t like = 0;
00145
00146 for (Int_t i=1; i<=fdPred->GetNbinsX(); ++i){
00147 //Bizarre limits because root histograms are silly
00148 Double_t mnu = fdPred->GetBinContent(i);
00149 Double_t dnu = fdData->GetBinContent(i);
00150 // if (i<10){
00151 // cout << "Bin " << i
00152 // << ": MC " << mnu
00153 // << "; data " << dnu
00154 // << endl;
00155 // }
00156 if (mnu<0.0001){mnu = 0.0001;}
00157 if (dnu){
00158 like += 2*(mnu - dnu + dnu*log(dnu/mnu));
00159 }
00160 else{like += 2*(mnu - dnu);}
00161 }
00162 // cout << "Like = " << like << endl;
00163 return like;
00164 }
|
|
|
Definition at line 292 of file NuMatrixFitterMinuit.cxx. References DoCPTFit(), ffdDataFilename, ffluxPoT, fhelperFilename, fndDataFilename, fxsecFilename, NuMatrixMethod::GetFDOtherCCTrueUnoscBackground(), NuMatrixMethod::PredictFDFluxUnosc(), NuMatrixMethod::SetExtrapolationScheme(), NuMatrixMethod::SetFDCCTrueUnoscBackground(), and NuMatrixMethod::WriteInputForFitter(). 00293 {
00294 //Get the FD data & PoT histograms
00295 //Data:
00296 TFile fdDataFile(ffdDataFilename.c_str(),"READ");
00297 TH1D* numufdData = (TH1D*) fdDataFile.Get("RecoEnergy_FD");
00298 numufdData->SetDirectory(0);
00299 numufdData->Sumw2();
00300 TH1D* numubarfdData = (TH1D*) fdDataFile.Get("RecoEnergyPQ_FD");
00301 numubarfdData->SetDirectory(0);
00302 numubarfdData->Sumw2();
00303
00304 //PoT:
00305 TH1F* hfdPoTTorTGT = (TH1F*) fdDataFile.Get("hPottortgt");
00306 hfdPoTTorTGT->SetDirectory(0);
00307 TH1F* hfdSpillsPerFile = (TH1F*) fdDataFile.Get("hSpillsPerFile");
00308 hfdSpillsPerFile->SetDirectory(0);
00309 TH1F* hfdRun = (TH1F*) fdDataFile.Get("hRun");
00310 hfdRun->SetDirectory(0);
00311 fdDataFile.Close();
00312
00313 //Get the ND PoT histograms
00314 TFile ndDataFile(fndDataFilename.c_str(),"READ");
00315 TH1F* hndPoTTorTGT = (TH1F*) ndDataFile.Get("hPottortgt");
00316 hndPoTTorTGT->SetDirectory(0);
00317 ndDataFile.Close();
00318
00319 //Calculate the total PoT
00320 Double_t ndPoT = hndPoTTorTGT->GetMean()*hndPoTTorTGT->GetEntries()*1e12;
00321 Double_t fdPoT = hfdRun->GetEntries()*hfdPoTTorTGT->GetMean()*1e8*1e12;
00322 // Double_t fdPoT = 6.5e20;
00323 cout << "ndPoT: " << ndPoT
00324 << "; fdPoT: " << fdPoT
00325 << endl;
00326
00327 //Hard coded, unchecked fiducial masses. Should make configurable.
00328 Double_t ndFidMass = 45.128*Munits::tonne;
00329 Double_t fdFidMass = 5.487*Munits::kilotonne*14.0*TMath::Pi()*448/(52.92*486);
00330
00331 NuMatrixMethod *fmmNuMu = new
00332 NuMatrixMethod(400,
00333 0,
00334 200,
00335 ndPoT,
00336 fdPoT,
00337 ffluxPoT,
00338 ndFidMass,
00339 fdFidMass,
00340 fhelperFilename,
00341 fxsecFilename,
00342 fndDataFilename,
00343 "$SRT_PRIVATE_CONTEXT/NtupleUtils/results/MMOutputMinuitNuMu.root");
00344 fmmNuMu->SetExtrapolationScheme(NuMMExtrapolation::kCCNuMuOscBackgroundCombined);
00345
00346 NuMatrixMethod *fmmNuMuBar = new
00347 NuMatrixMethod(400,
00348 0,
00349 200,
00350 ndPoT,
00351 fdPoT,
00352 ffluxPoT,
00353 ndFidMass,
00354 fdFidMass,
00355 fhelperFilename,
00356 fxsecFilename,
00357 fndDataFilename,
00358 "$SRT_PRIVATE_CONTEXT/NtupleUtils/results/MMOutputMinuitNuMuBar.root");
00359 fmmNuMuBar->SetExtrapolationScheme(NuMMExtrapolation::kCCNuMuBarOscBackgroundCombined);
00360
00361 fmmNuMu->PredictFDFluxUnosc();
00362 fmmNuMuBar->PredictFDFluxUnosc();
00363
00364 fmmNuMu->SetFDCCTrueUnoscBackground(fmmNuMuBar->GetFDOtherCCTrueUnoscBackground());
00365 fmmNuMuBar->SetFDCCTrueUnoscBackground(fmmNuMu->GetFDOtherCCTrueUnoscBackground());
00366
00367 NuMatrixInput *mmInput = new NuMatrixInput();
00368
00369 fmmNuMu->WriteInputForFitter(mmInput);
00370 fmmNuMuBar->WriteInputForFitter(mmInput);
00371
00372 /*NuMatrixOutput *mmOut = */this->DoCPTFit(mmInput, numufdData, numubarfdData);
00373
00374 delete fmmNuMu; fmmNuMu = 0;
00375 delete fmmNuMuBar; fmmNuMuBar = 0;
00376 }
|
|
||||||||||||||||
|
Definition at line 168 of file NuMatrixFitterMinuit.cxx. References DoCPTFit(), NuXMLConfig::FullTitle(), infile, and NuMatrixOutput::XMLConfig(). 00173 {
00174 TFile fdDataFile(fdDataFilename.c_str(),"READ");
00175 TH1D* numubarfdData = (TH1D*) fdDataFile.Get("RecoEnergyPQ_FD");
00176 numubarfdData->SetDirectory(0);
00177 numubarfdData->Sumw2();
00178
00179 TH1D* numufdData = (TH1D*) fdDataFile.Get("RecoEnergy_FD");
00180 numufdData->SetDirectory(0);
00181 numufdData->Sumw2();
00182 fdDataFile.Close();
00183
00184
00185 TFile *infile = new TFile(inputFilename.c_str(), "READ");
00186
00187
00188 NuMatrixInput *mmInput = 0;
00189 mmInput = (NuMatrixInput*) infile->Get("NuMatrixInput");
00190 if (!mmInput) {
00191 mmInput = new NuMatrixInput(infile);
00192 }
00193 NuMatrixOutput *mmOutput = this->DoCPTFit(mmInput, numufdData, numubarfdData);
00194
00195 TString objName = "NuMatrixOutput";
00196 NuXMLConfig *xmlConfig = (NuXMLConfig*)gROOT->FindObject("NuXMLConfig");
00197
00198 if (xmlConfig) {
00199 cout << "This file fit is for systematic: " << xmlConfig->FullTitle() << endl;
00200 mmOutput->XMLConfig(xmlConfig);
00201 //objName = xmlConfig->FullName();
00202 }
00203 else {
00204 cout << "No NuXMLConfig object -- nominal fit assumed." << endl;
00205 //objName = "Nominal";
00206 }
00207 cout << "Recreating file " << outputFilename << endl;
00208 TFile *outputFile = new TFile(outputFilename.c_str(),"RECREATE");
00209 outputFile->cd();
00210 mmOutput->Write(objName);
00211 outputFile->Close();
00212 infile->cd();
00213 cout << "Don't mind the segfault. Everything is done and it appears to be harmless." << endl;
00214 }
|
|
||||||||||||||||
|
Definition at line 218 of file NuMatrixFitterMinuit.cxx. References NuMatrixOutput::BestFitPoint(), fFitter, fInput, fnumubarfdData, fnumufdData, NuMatrixOutput::NuMuBarBestFitFDPrediction(), NuMatrixOutput::NuMuBarFDData(), NuMatrixOutput::NuMuBarNoOscPrediction(), NuMatrixInput::PredictFDSpectrumNuMuBar(), and SetupPars(). Referenced by CPTAnalysis(), and CPTFit(). 00219 {
00220 fnumufdData = numufdData;
00221 fnumubarfdData = numubarfdData;
00222
00223 fnumufdData->Sumw2();
00224 fnumubarfdData->Sumw2();
00225
00226 //fmmNuMu = new NuMatrixMethod(mmInput,NuMMExtrapolation::kModularNuMuCPTFit);
00227 //fmmNuMuBar = new NuMatrixMethod(mmInput,NuMMExtrapolation::kModularNuMuBarCPTFit);
00228 fInput = mmInput;
00229
00230 this->SetupPars();
00231 fFitter->Minimize();
00232 fFitter->PrintResults(0,0);
00233
00234 Double_t dm2 = fFitter->GetParameter(0);
00235 Double_t sn2 = fFitter->GetParameter(1);
00236 Double_t dm2bar = fFitter->GetParameter(2);
00237 Double_t sn2bar = fFitter->GetParameter(3);
00238
00239 NuMatrixOutput *mmOutput = new NuMatrixOutput();
00240 //mmOutput->NuMuBarChi2Surface(hChi2Bar);
00241 mmOutput->NuMuBarFDData(numubarfdData);
00242
00243 mmOutput->NuMuBarBestFitFDPrediction(*fInput->PredictFDSpectrumNuMuBar(dm2, sn2, dm2bar, sn2bar));
00244 mmOutput->NuMuBarNoOscPrediction(*fInput->PredictFDSpectrumNuMuBar(0., 0., 0., 0.));
00245
00246 mmOutput->BestFitPoint(sn2bar, dm2bar);
00247
00248 return mmOutput;
00249 //delete fmmNuMu; fmmNuMu = 0;
00250 //delete fmmNuMuBar; fmmNuMuBar = 0;
00251 }
|
|
|
Definition at line 254 of file NuMatrixFitterMinuit.cxx. References CalculateLikelihood(), fInput, fnumubarfdData, fnumufdData, NuMatrixInput::PredictFDSpectrumNuMu(), and NuMatrixInput::PredictFDSpectrumNuMuBar(). 00255 {
00256 static Int_t fitcounter = 0;
00257 if (!(fitcounter%50)){
00258 cout << "Minuit TF call " << fitcounter << endl;
00259 }
00260 ++fitcounter;
00261
00262 Double_t dm2 = pars[0];
00263 Double_t sn2 = pars[1];
00264 Double_t dm2bar = pars[2];
00265 Double_t sn2bar = pars[3];
00266 // Double_t normalisation = pars[4];
00267
00268 TH1D* numufdPred = fInput->PredictFDSpectrumNuMu(dm2,sn2,dm2bar,sn2bar);
00269 TH1D* numubarfdPred = fInput->PredictFDSpectrumNuMuBar(dm2,sn2,dm2bar,sn2bar);
00270
00271 // numufdPred->Scale(normalisation);
00272 // numubarfdPred->Scale(normalisation);
00273 // cout << "NuMu fit:" << endl;
00274 Double_t chi2 = this->CalculateLikelihood(numufdPred,fnumufdData);
00275 // cout << "NuBar fit:" << endl;
00276 chi2 += this->CalculateLikelihood(numubarfdPred,fnumubarfdData);
00277
00278 // cout << "dm2: " << dm2 << "; sn2: " << sn2
00279 // << "; dm2bar: " << dm2bar << "; sn2bar: " << sn2bar
00280 // << "; chi2: " << chi2 << endl;
00281
00282 return chi2;
00283 }
|
|
|
Definition at line 94 of file NuMatrixFitterMinuit.cxx. References fFitter. Referenced by DoCPTFit(). 00095 {
00096 fFitter->SetParameter(0,
00097 "Dm2",
00098 0.003,0.2,1.0,0.0);
00099 fFitter->SetParameter(1,
00100 "Sn2",
00101 1.0,0.2,1.0,0.0);
00102 fFitter->SetParameter(2,
00103 "Dm2Bar",
00104 0.003,0.2,1.0,0.0);
00105 fFitter->SetParameter(3,
00106 "Sn2Bar",
00107 1.0,0.2,1.0,0.0);
00108 fFitter->SetParameter(4,
00109 "Normalisation",
00110 1.0,0.2,1.0,0.0);
00111 fFitter->FixParameter(4);
00112 }
|
|
||||||||||||||||||||
|
Definition at line 379 of file NuMatrixFitterMinuit.cxx. References CalculateLikelihood(), ffdDataFilename, ffluxPoT, fhelperFilename, fndDataFilename, fnumubarfdData, fnumufdData, fxsecFilename, NuMatrixMethod::GetFDOtherCCTrueUnoscBackground(), NuMatrixMethod::PredictFDFluxUnosc(), NuMatrixMethod::PredictFDSpectrumBackgroundOscCombined(), NuMatrixMethod::SetExtrapolationScheme(), and NuMatrixMethod::SetFDCCTrueUnoscBackground(). 00383 {
00384 //Get the FD data & PoT histograms
00385 //Data:
00386 TFile fdDataFile(ffdDataFilename.c_str(),"READ");
00387 fnumufdData = (TH1D*) fdDataFile.Get("RecoEnergy_FD");
00388 fnumufdData->SetDirectory(0);
00389 fnumufdData->Sumw2();
00390 fnumubarfdData = (TH1D*) fdDataFile.Get("RecoEnergyPQ_FD");
00391 fnumubarfdData->SetDirectory(0);
00392 fnumubarfdData->Sumw2();
00393 //PoT:
00394 TH1F* hfdPoTTorTGT = (TH1F*) fdDataFile.Get("hPottortgt");
00395 hfdPoTTorTGT->SetDirectory(0);
00396 TH1F* hfdSpillsPerFile = (TH1F*) fdDataFile.Get("hSpillsPerFile");
00397 hfdSpillsPerFile->SetDirectory(0);
00398 TH1F* hfdRun = (TH1F*) fdDataFile.Get("hRun");
00399 hfdRun->SetDirectory(0);
00400 fdDataFile.Close();
00401
00402 //Get the ND PoT histograms
00403 TFile ndDataFile(fndDataFilename.c_str(),"READ");
00404 TH1F* hndPoTTorTGT = (TH1F*) ndDataFile.Get("hPottortgt");
00405 hndPoTTorTGT->SetDirectory(0);
00406 ndDataFile.Close();
00407
00408 //Calculate the total PoT
00409 Double_t ndPoT = hndPoTTorTGT->GetMean()*hndPoTTorTGT->GetEntries()*1e12;
00410 Double_t fdPoT = hfdRun->GetEntries()*hfdPoTTorTGT->GetMean()*1e8*1e12;
00411 // Double_t fdPoT = 6.5e20;
00412 cout << "ndPoT: " << ndPoT
00413 << "; fdPoT: " << fdPoT
00414 << endl;
00415
00416 //Hard coded, unchecked fiducial masses. Should make configurable.
00417 Double_t ndFidMass = 45.128*Munits::tonne;
00418 Double_t fdFidMass = 5.487*Munits::kilotonne*14.0*TMath::Pi()*448/(52.92*486);
00419
00420 NuMatrixMethod *fmmNuMu = new
00421 NuMatrixMethod(400,
00422 0,
00423 200,
00424 ndPoT,
00425 fdPoT,
00426 ffluxPoT,
00427 ndFidMass,
00428 fdFidMass,
00429 fhelperFilename,
00430 fxsecFilename,
00431 fndDataFilename,
00432 "$SRT_PRIVATE_CONTEXT/NtupleUtils/results/MMOutputMinuitNuMu.root");
00433 fmmNuMu->SetExtrapolationScheme
00434 (NuMMExtrapolation::kCCNuMuOscBackgroundCombined);
00435
00436 NuMatrixMethod *fmmNuMuBar = new
00437 NuMatrixMethod(400,
00438 0,
00439 200,
00440 ndPoT,
00441 fdPoT,
00442 ffluxPoT,
00443 ndFidMass,
00444 fdFidMass,
00445 fhelperFilename,
00446 fxsecFilename,
00447 fndDataFilename,
00448 "$SRT_PRIVATE_CONTEXT/NtupleUtils/results/MMOutputMinuitNuMuBar.root");
00449 fmmNuMuBar->SetExtrapolationScheme
00450 (NuMMExtrapolation::kCCNuMuBarOscBackgroundCombined);
00451
00452 fmmNuMu->PredictFDFluxUnosc();
00453 fmmNuMuBar->PredictFDFluxUnosc();
00454
00455 fmmNuMu->SetFDCCTrueUnoscBackground(fmmNuMuBar->GetFDOtherCCTrueUnoscBackground());
00456 fmmNuMuBar->SetFDCCTrueUnoscBackground(fmmNuMu->GetFDOtherCCTrueUnoscBackground());
00457
00458
00459 const TH1D* numufdPred =
00460 fmmNuMu->PredictFDSpectrumBackgroundOscCombined(dm2,
00461 sn2,
00462 dm2bar,
00463 sn2bar);
00464 const TH1D* numubarfdPred =
00465 fmmNuMuBar->PredictFDSpectrumBackgroundOscCombined(dm2bar,
00466 sn2bar,
00467 dm2,
00468 sn2);
00469 // cout << "NuMu fit:" << endl;
00470 Double_t chi2 = this->CalculateLikelihood(numufdPred,fnumufdData);
00471 // cout << "NuBar fit:" << endl;
00472 chi2 += this->CalculateLikelihood(numubarfdPred,fnumubarfdData);
00473
00474 // cout << "dm2: " << dm2 << "; sn2: " << sn2
00475 // << "; dm2bar: " << dm2bar << "; sn2bar: " << sn2bar
00476 // << "; chi2: " << chi2 << endl;
00477
00478 delete fmmNuMu; fmmNuMu = 0;
00479 delete fmmNuMuBar; fmmNuMuBar = 0;
00480 TFile NuMuFile("$SRT_PRIVATE_CONTEXT/NtupleUtils/results/MMOutputMinuitNuMu.root","UPDATE");
00481 fnumufdData->Write();
00482 NuMuFile.Close();
00483 TFile NuMuBarFile("$SRT_PRIVATE_CONTEXT/NtupleUtils/results/MMOutputMinuitNuMuBar.root","UPDATE");
00484 fnumubarfdData->Write();
00485 NuMuBarFile.Close();
00486 }
|
|
|
Definition at line 45 of file NuMatrixFitterMinuit.h. 00045 {return 1.0;}//Chi2 errors; 0.5 for likelihood
|
|
|
Definition at line 71 of file NuMatrixFitterMinuit.h. Referenced by CPTAnalysis(), NuMatrixFitterMinuit(), and SingleCPTPrediction(). |
|
|
Definition at line 63 of file NuMatrixFitterMinuit.h. Referenced by DoCPTFit(), NuMatrixFitterMinuit(), and SetupPars(). |
|
|
Definition at line 69 of file NuMatrixFitterMinuit.h. Referenced by CPTAnalysis(), NuMatrixFitterMinuit(), and SingleCPTPrediction(). |
|
|
Definition at line 72 of file NuMatrixFitterMinuit.h. Referenced by CPTAnalysis(), NuMatrixFitterMinuit(), and SingleCPTPrediction(). |
|
|
Definition at line 62 of file NuMatrixFitterMinuit.h. Referenced by DoCPTFit(), NuMatrixFitterMinuit(), and operator()(). |
|
|
Definition at line 70 of file NuMatrixFitterMinuit.h. Referenced by CPTAnalysis(), NuMatrixFitterMinuit(), and SingleCPTPrediction(). |
|
|
Definition at line 65 of file NuMatrixFitterMinuit.h. Referenced by DoCPTFit(), operator()(), and SingleCPTPrediction(). |
|
|
Definition at line 64 of file NuMatrixFitterMinuit.h. Referenced by DoCPTFit(), operator()(), and SingleCPTPrediction(). |
|
|
Definition at line 74 of file NuMatrixFitterMinuit.h. Referenced by NuMatrixFitterMinuit(). |
|
|
Definition at line 73 of file NuMatrixFitterMinuit.h. Referenced by CPTAnalysis(), NuMatrixFitterMinuit(), and SingleCPTPrediction(). |
1.3.9.1