loon: NuHistos.cxx Source File

00001 
00002 
00003 // Coded by Jeff Hartnell Jul/2005
00004 //
00005 // Contact: j.j.hartnell@rl.ac.uk
00007 
00008 #include <cassert>
00009 #include <cmath>
00010 
00011 #include "TH2D.h"
00012 #include "TMath.h"
00013 #include "TROOT.h"
00014 
00015 #include "Conventions/Detector.h"
00016 #include "MessageService/MsgService.h"
00017 #include "Conventions/SimFlag.h"
00018 
00019 #include "NtupleUtils/NuEvent.h"
00020 #include "NtupleUtils/NuCuts.h"
00021 #include "NtupleUtils/NuHistos.h"
00022 
00023 using std::endl;
00024 using std::cout;
00025 
00026 using namespace ENCTruth;
00027 
00028 CVSID("$Id: NuHistos.cxx,v 1.38 2009/10/13 17:29:11 lefeuvre Exp $");
00029 
00030 //......................................................................
00031 
00032 NuHistos::NuHistos()
00033 {
00034   MSG("NuHistos",Msg::kDebug)
00035     <<"Running NuHistos Constructor..."<<endl;
00036 
00037 
00038   MSG("NuHistos",Msg::kDebug)
00039     <<"Finished NuHistos Constructor"<<endl;
00040 }
00041 
00042 //......................................................................
00043 
00044 NuHistos::~NuHistos()
00045 {
00046   MSG("NuHistos",Msg::kDebug)
00047     <<"Running NuHistos Destructor..."<<endl;
00048 
00049 
00050   MSG("NuHistos",Msg::kDebug)
00051     <<"Finished NuHistos Destructor"<<endl;
00052 }
00053 
00054 //......................................................................
00055 
00056 void NuHistos::FillMatrixMethodHistos
00057 (const NuEvent& nu,
00058  const NuBinningScheme::NuBinningScheme_t binningScheme) const
00059 {
00060   //TH1D *fFracErrOnPred;//JE does this one
00061   //TH1D *fTrueEnergyTrueCCFidEvents_ND;//this is in NuPlots
00062   //TH1D *fTrueEnergyTrueCCFidEvents_FD;//this is in NuPlots
00063 
00064   MAXMSG("NuHistos",Msg::kInfo,1)
00065     <<"Running FillMatrixMethodHistos..."<<endl;
00066 
00067   //MC Helper histos (CC numus):
00068   static TH2D* fRecoVsTrueEnergy_ND=0;
00069   static TH1D* fEfficiency_ND=0;
00070   static TH1D* fPurity_ND=0;
00071   static TH2D* fRecoVsTrueEnergy_FD=0;
00072   static TH1D* fEfficiency_FD=0;
00073   static TH1D* fPurity_FD=0;
00074   static TH1D* fCCContamination_ND = 0;
00075   static TH1D* fCCContamination_FD = 0;
00076   static TH2D* fCCContaminationRecoVsTrue_FD = 0;
00077   static TH1D* fNCContamination_ND = 0;
00078   static TH1D* fNCContamination_FD = 0;
00079   //MC Helper histos (CC numubars):
00080   static TH2D* fRecoVsTrueEnergyPQ_ND=0;
00081   static TH1D* fEfficiencyPQ_ND=0;
00082   static TH1D* fPurityPQ_ND=0;
00083   static TH2D* fRecoVsTrueEnergyPQ_FD=0;
00084   static TH1D* fEfficiencyPQ_FD=0;
00085   static TH1D* fPurityPQ_FD=0;
00086   static TH1D* fCCContaminationPQ_ND = 0;
00087   static TH1D* fCCContaminationPQ_FD = 0;
00088   static TH2D* fCCContaminationRecoVsTruePQ_FD = 0;
00089   static TH1D* fNCContaminationPQ_ND = 0;
00090   static TH1D* fNCContaminationPQ_FD = 0;
00091   //MC Helper histos (CC numus + CC numubars):
00092   static TH2D* fRecoVsTrueEnergyAll_ND=0;
00093   static TH1D* fEfficiencyAll_ND=0;
00094   static TH1D* fPurityAll_ND=0;
00095   static TH2D* fRecoVsTrueEnergyAll_FD=0;
00096   static TH1D* fEfficiencyAll_FD=0;
00097   static TH1D* fPurityAll_FD=0;
00098   //MC Helper histos (CC nutaus):
00099   static TH2D* fRecoVsTrueEnergyTau_FD=0;
00100   static TH1D* fEfficiencyTau_FD=0;
00101   //MC Helper histos (CC nutaus):
00102   static TH2D* fRecoVsTrueEnergyTauPQ_FD=0;
00103   static TH1D* fEfficiencyTauPQ_FD=0;
00104   //MC Helper histos (CC nutaus):
00105   static TH2D* fRecoVsTrueEnergyTauAll_FD=0;
00106   static TH1D* fEfficiencyTauAll_FD=0;
00107 
00108   //MC Check list histos (CC numus):
00109   static TH1D* fTrueEnergyCCOnlyEvents_ND=0;
00110   static TH1D* fRecoEnergyCCOnlyEvents_ND=0;
00111   static TH1D* fRecoEnergyAllEvents_ND=0;
00112   static TH1D* fTrueEnergyCCOnlyEvents_FD=0;
00113   static TH1D* fRecoEnergyCCOnlyEvents_FD=0;
00114   static TH1D* fRecoEnergyAllEvents_FD=0;
00115   //MC Check list histos (CC numubars):
00116   static TH1D* fTrueEnergyCCOnlyEventsPQ_ND=0;
00117   static TH1D* fRecoEnergyCCOnlyEventsPQ_ND=0;
00118   static TH1D* fRecoEnergyAllEventsPQ_ND=0;
00119   static TH1D* fTrueEnergyCCOnlyEventsPQ_FD=0;
00120   static TH1D* fRecoEnergyCCOnlyEventsPQ_FD=0;
00121   static TH1D* fRecoEnergyAllEventsPQ_FD=0;
00122   //MC Check list histos (CC numus + CC numubars):
00123   static TH1D* fTrueEnergyCCOnlyEventsAll_ND=0;
00124   static TH1D* fRecoEnergyCCOnlyEventsAll_ND=0;
00125   static TH1D* fRecoEnergyAllEventsAll_ND=0;
00126   static TH1D* fTrueEnergyCCOnlyEventsAll_FD=0;
00127   static TH1D* fRecoEnergyCCOnlyEventsAll_FD=0;
00128   static TH1D* fRecoEnergyAllEventsAll_FD=0;
00129   //MC checklist histos (CC nutaus):
00130   static TH1D* fTrueEnergyCCOnlyEventsTau_FD=0;
00131   //MC checklist histos (CC nutaubars):
00132   static TH1D* fTrueEnergyCCOnlyEventsTauPQ_FD=0;
00133   //MC checklist histos (CC nutaus + nutaubars):
00134   static TH1D* fTrueEnergyCCOnlyEventsTauAll_FD=0;
00135 
00136   //DATA histos (CC numus)
00137   static TH1D* RecoEnergy_ND=0;
00138   static TH1D* RecoEnergy_FD=0;
00139   //DATA histos (CC numubars)
00140   static TH1D* RecoEnergyPQ_ND=0;
00141   static TH1D* RecoEnergyPQ_FD=0;
00142 
00143   //DATA histos (Pure Reco Energys)
00144   static TH1D* RecoEnergyTrueTauPQ_FD=0;
00145   static TH1D* RecoEnergyTrueMuonPQ_FD=0;
00146   static TH1D* RecoEnergyTrueTauNQ_FD=0;
00147   static TH1D* RecoEnergyTrueMuonNQ_FD=0;
00148 
00149   static TH1D* RecoEnergyTrueTauPQ_ND=0;
00150   static TH1D* RecoEnergyTrueMuonPQ_ND=0;
00151   static TH1D* RecoEnergyTrueTauNQ_ND=0;
00152   static TH1D* RecoEnergyTrueMuonNQ_ND=0;
00153 
00154   static TH1D* RecoEnergyWSTauPQ_FD=0;
00155   static TH1D* RecoEnergyWSMuonPQ_FD=0;
00156   static TH1D* RecoEnergyWSTauNQ_FD=0;
00157   static TH1D* RecoEnergyWSMuonNQ_FD=0;
00158 
00159   static TH1D* RecoEnergyWSTauPQ_ND=0;
00160   static TH1D* RecoEnergyWSMuonPQ_ND=0;
00161   static TH1D* RecoEnergyWSTauNQ_ND=0;
00162   static TH1D* RecoEnergyWSMuonNQ_ND=0;
00163 
00164 
00165   //DATA histos (CC numus + CC numubars)
00166   static TH1D* RecoEnergyAll_ND=0;
00167   static TH1D* RecoEnergyAll_FD=0;
00168 
00169   static TH1F*  poKinPQ = 0;
00170   static TH1F*  roIdPQ = 0;
00171   static TH1F*  dpIdPQ = 0;
00172   static TH1F*  trkLPQ = 0;
00173   static TH1F*  relAPQ = 0;
00174   static TH1F*  relANQ = 0;
00175   static TH1F*  relAPQ1 = 0;
00176   static TH1F*  qp_sigmaqpPQ = 0;
00177   static TH1F*  majCPQ = 0;
00178   static TH1F*  probPQ = 0;
00179   static TH2F*  relAngVsTrkEndX = 0;
00180   static TH2F*  relAngVsTrkEndY = 0;
00181   static TH2F*  relAngVsTrkEndZ = 0;
00182 
00183   if(!fRecoVsTrueEnergy_ND){
00184     MSG("NuHistos",Msg::kInfo)
00185       <<"Creating histograms on first run of FillMatrixMethodHistos..."
00186       <<" gDirectory gives:"<<endl;
00187           gDirectory->Print();
00188 
00189     const NuUtilities cuts;
00190     const vector<Double_t> recoBins = cuts.RecoBins(binningScheme);
00191     const Int_t numRecoBins = recoBins.size()-1;
00192     Float_t *recoBinsArray;
00193     recoBinsArray=new Float_t[numRecoBins+1];
00194     {
00195       Int_t i=0;
00196       for (vector<Double_t>::const_iterator itBin = recoBins.begin();
00197            itBin != recoBins.end();
00198            ++itBin, ++i){
00199         recoBinsArray[i] = *itBin;
00200       }
00201     }
00202 
00203     const vector<Double_t> trueBins = cuts.TrueBins(binningScheme);
00204     const Int_t numTrueBins = trueBins.size()-1;
00205     Float_t *trueBinsArray;
00206     trueBinsArray=new Float_t[numTrueBins+1];
00207     {
00208       Int_t i=0;
00209       for (vector<Double_t>::const_iterator itBin = trueBins.begin();
00210            itBin != trueBins.end();
00211            ++itBin, ++i){
00212         trueBinsArray[i] = *itBin;
00213       }
00214     }
00215 
00216     //helper histos (CC numus):
00217     fRecoVsTrueEnergy_ND=new TH2D
00218                 ("RecoVsTrueEnergy_ND",
00219                  "Reco Vs True Energy (NearDet)",
00220                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00221     fRecoVsTrueEnergy_ND->Sumw2();
00222 
00223     fEfficiency_ND=new TH1D
00224                 ("Efficiency_ND",
00225                  "NuMu CC Selection Efficiency with True Energy (NearDet)",
00226                  numTrueBins,trueBinsArray);
00227     fEfficiency_ND->Sumw2();
00228 
00229     fPurity_ND=new TH1D
00230                 ("Purity_ND",
00231                  "NuMu CC Selection Purity with Reco Energy (NearDet)",
00232                  numRecoBins,recoBinsArray);
00233     fPurity_ND->Sumw2();
00234 
00235     fRecoVsTrueEnergy_FD=new TH2D
00236                 ("RecoVsTrueEnergy_FD",
00237                  "Reco Vs True Energy (FarDet)",
00238                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00239     fRecoVsTrueEnergy_FD->Sumw2();
00240 
00241     fEfficiency_FD=new TH1D
00242                 ("Efficiency_FD",
00243                  "NuMu CC Selection Efficiency with True Energy (FarDet)",
00244                  numTrueBins,trueBinsArray);
00245     fEfficiency_FD->Sumw2();
00246 
00247     fPurity_FD=new TH1D
00248                 ("Purity_FD",
00249                  "NuMu CC Selection Purity with Reco Energy (FarDet)",
00250                  numRecoBins,recoBinsArray);
00251     fPurity_FD->Sumw2();
00252 
00253     fCCContamination_FD = new TH1D
00254                 ("CCContamination_FD",
00255                  "CC antineutrino contamination in a -ve curving CC selection (FarDet)",
00256                  numTrueBins,trueBinsArray);
00257     fCCContamination_FD->Sumw2();
00258 
00259     fCCContamination_ND = new TH1D
00260                 ("CCContamination_ND",
00261                  "CC antineutrino contamination in a -ve curving CC selection (NearDet)",
00262                  numRecoBins,recoBinsArray);
00263     fCCContamination_ND->Sumw2();
00264 
00265     fCCContaminationRecoVsTrue_FD = new TH2D
00266                 ("CCContaminationRecoVsTrue_FD",
00267                  "Reco vs true energy of CC antineutrino contamination in a -ve curving CC selection (FarDet)",
00268                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00269     fCCContaminationRecoVsTrue_FD->Sumw2();
00270 
00271     fNCContamination_ND = new TH1D
00272                 ("NCContamination_ND",
00273                  "NC contamination in a -ve curving CC selection (NearDet)",
00274                  numRecoBins,recoBinsArray);
00275     fNCContamination_ND->Sumw2();
00276 
00277     fNCContamination_FD = new TH1D
00278                 ("NCContamination_FD",
00279                  "NC contamination in a -ve curving CC selection (FarDet)",
00280                  numRecoBins,recoBinsArray);
00281     fNCContamination_FD->Sumw2();
00282 
00283     //helper histos (CC numubars):
00284     fRecoVsTrueEnergyPQ_ND=new TH2D
00285                 ("RecoVsTrueEnergyPQ_ND",
00286                  "Reco Vs True Energy (NearDet)",
00287                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00288     fRecoVsTrueEnergyPQ_ND->Sumw2();
00289 
00290     fEfficiencyPQ_ND=new TH1D
00291                 ("EfficiencyPQ_ND",
00292                  "NuMuBar CC Selection Efficiency with True Energy (NearDet)",
00293                  numTrueBins,trueBinsArray);
00294     fEfficiencyPQ_ND->Sumw2();
00295 
00296     fPurityPQ_ND=new TH1D
00297                 ("PurityPQ_ND",
00298                  "NuMuBar CC Selection Purity with Reco Energy (NearDet)",
00299                  numRecoBins,recoBinsArray);
00300     fPurityPQ_ND->Sumw2();
00301 
00302     fRecoVsTrueEnergyPQ_FD=new TH2D
00303                 ("RecoVsTrueEnergyPQ_FD",
00304                  "Reco Vs True Energy (FarDet)",
00305                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00306     fRecoVsTrueEnergyPQ_FD->Sumw2();
00307 
00308     fEfficiencyPQ_FD=new TH1D
00309                 ("EfficiencyPQ_FD",
00310                  "NuMuBar CC Selection Efficiency with True Energy (FarDet)",
00311                  numTrueBins,trueBinsArray);
00312     fEfficiencyPQ_FD->Sumw2();
00313 
00314     fPurityPQ_FD=new TH1D
00315                 ("PurityPQ_FD",
00316                  "NuMuBar CC Selection Purity with Reco Energy (FarDet)",
00317                  numRecoBins,recoBinsArray);
00318     fPurityPQ_FD->Sumw2();
00319 
00320     fCCContaminationPQ_ND = new TH1D
00321                 ("CCContaminationPQ_ND",
00322                  "CC neutrino contamination in a +ve curving CC selection (NearDet)",
00323                  numRecoBins,recoBinsArray);
00324     fCCContaminationPQ_ND->Sumw2();
00325 
00326     fCCContaminationPQ_FD = new TH1D
00327                 ("CCContaminationPQ_FD",
00328                  "CC neutrino contamination in a +ve curving CC selection (FarDet)",
00329                  numTrueBins,trueBinsArray);
00330     fCCContaminationPQ_FD->Sumw2();
00331 
00332     fCCContaminationRecoVsTruePQ_FD = new TH2D
00333                 ("CCContaminationRecoVsTruePQ_FD",
00334                  "RecoVsTrue Energy of CC neutrino contamination in a +ve curving CC selection (FarDet)",
00335                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00336     fCCContaminationRecoVsTruePQ_FD->Sumw2();
00337 
00338     fNCContaminationPQ_ND = new TH1D
00339                 ("NCContaminationPQ_ND",
00340                  "NC contamination in a +ve curving CC selection (NearDet)",
00341                  numRecoBins,recoBinsArray);
00342     fNCContaminationPQ_ND->Sumw2();
00343 
00344     fNCContaminationPQ_FD = new TH1D
00345                 ("NCContaminationPQ_FD",
00346                  "NC contamination in a +ve curving CC selection (FarDet)",
00347                  numRecoBins,recoBinsArray);
00348     fNCContaminationPQ_FD->Sumw2();
00349 
00350     //helper histos (CC numus + CC numubars):
00351     fRecoVsTrueEnergyAll_ND=new TH2D
00352                 ("RecoVsTrueEnergyAll_ND",
00353                  "Reco Vs True Energy (NearDet)",
00354                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00355     fRecoVsTrueEnergyAll_ND->Sumw2();
00356 
00357     fEfficiencyAll_ND=new TH1D
00358                 ("EfficiencyAll_ND",
00359                  "All CC Selection Efficiency with True Energy (NearDet)",
00360                  numTrueBins,trueBinsArray);
00361     fEfficiencyAll_ND->Sumw2();
00362 
00363     fPurityAll_ND=new TH1D
00364                 ("PurityAll_ND",
00365                  "All CC Selection Purity with Reco Energy (NearDet)",
00366                  numRecoBins,recoBinsArray);
00367     fPurityAll_ND->Sumw2();
00368 
00369     fRecoVsTrueEnergyAll_FD=new TH2D
00370                 ("RecoVsTrueEnergyAll_FD",
00371                  "Reco Vs True Energy (FarDet)",
00372                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00373     fRecoVsTrueEnergyAll_FD->Sumw2();
00374 
00375     fEfficiencyAll_FD=new TH1D
00376                 ("EfficiencyAll_FD",
00377                  "All CC Selection Efficiency with True Energy (FarDet)",
00378                  numTrueBins,trueBinsArray);
00379     fEfficiencyAll_FD->Sumw2();
00380 
00381     fPurityAll_FD=new TH1D
00382                 ("PurityAll_FD",
00383                  "All CC Selection Purity with Reco Energy (FarDet)",
00384                  numRecoBins,recoBinsArray);
00385     fPurityAll_FD->Sumw2();
00386 
00387     //Helper histos (CC nutau):
00388     fRecoVsTrueEnergyTau_FD=new TH2D
00389                 ("RecoVsTrueEnergyTau_FD",
00390                  "Tau Reco Vs True Energy (FarDet)",
00391                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00392     fRecoVsTrueEnergyTau_FD->Sumw2();
00393 
00394     fEfficiencyTau_FD=new TH1D
00395                 ("EfficiencyTau_FD",
00396                  "NuTau CC Selection Efficiency with True Energy (FarDet)",
00397                  numTrueBins,trueBinsArray);
00398     fEfficiencyTau_FD->Sumw2();
00399 
00400     //Helper histos (CC nutaubar):
00401     fRecoVsTrueEnergyTauPQ_FD=new TH2D
00402                 ("RecoVsTrueEnergyTauPQ_FD",
00403                  "TauBar Reco Vs True Energy (FarDet)",
00404                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00405     fRecoVsTrueEnergyTauPQ_FD->Sumw2();
00406 
00407     fEfficiencyTauPQ_FD=new TH1D
00408                 ("EfficiencyTauPQ_FD",
00409                  "NuTauBar CC Selection Efficiency with True Energy (FarDet)",
00410                  numTrueBins,trueBinsArray);
00411     fEfficiencyTauPQ_FD->Sumw2();
00412 
00413     //Helper histos (CC nutau + nutaubar):
00414     fRecoVsTrueEnergyTauAll_FD=new TH2D
00415                 ("RecoVsTrueEnergyTauAll_FD",
00416                  "Tau Reco Vs True Energy (FarDet)",
00417                  numTrueBins,trueBinsArray,numRecoBins,recoBinsArray);
00418     fRecoVsTrueEnergyTauAll_FD->Sumw2();
00419 
00420     fEfficiencyTauAll_FD=new TH1D
00421                 ("EfficiencyTauAll_FD",
00422                  "NuTau CC Selection Efficiency with True Energy (FarDet)",
00423                  numTrueBins,trueBinsArray);
00424     fEfficiencyTauAll_FD->Sumw2();
00425 
00426     //Checklist histos (CC numus):
00427     fTrueEnergyCCOnlyEvents_ND=new TH1D
00428                 ("TrueEnergyCCOnlyEvents_ND",
00429                  "NuMu CC Selected Events with True Energy (NearDet)",
00430                  numTrueBins,trueBinsArray);
00431     fTrueEnergyCCOnlyEvents_ND->Sumw2();
00432 
00433     fRecoEnergyCCOnlyEvents_ND=new TH1D
00434                 ("RecoEnergyCCOnlyEvents_ND",
00435                  "NuMu CC Selected Events with Reco Energy (NearDet)",
00436                  numRecoBins,recoBinsArray);
00437     fRecoEnergyCCOnlyEvents_ND->Sumw2();
00438 
00439     fRecoEnergyAllEvents_ND=new TH1D
00440                 ("RecoEnergyAllEvents_ND",
00441                  "Selected Events with Reco Energy (NearDet)",
00442                  numRecoBins,recoBinsArray);
00443     fRecoEnergyAllEvents_ND->Sumw2();
00444 
00445     fTrueEnergyCCOnlyEvents_FD=new TH1D
00446                 ("TrueEnergyCCOnlyEvents_FD",
00447                  "NuMu CC Selected Events with True Energy (FarDet)",
00448                  numTrueBins,trueBinsArray);
00449     fTrueEnergyCCOnlyEvents_FD->Sumw2();
00450 
00451     fRecoEnergyCCOnlyEvents_FD=new TH1D
00452                 ("RecoEnergyCCOnlyEvents_FD",
00453                  "NuMu CC Selected Events with Reco Energy (FarDet)",
00454                  numRecoBins,recoBinsArray);
00455     fRecoEnergyCCOnlyEvents_FD->Sumw2();
00456 
00457     fRecoEnergyAllEvents_FD=new TH1D
00458                 ("RecoEnergyAllEvents_FD",
00459                  "Selected Events with Reco Energy (FarDet)",
00460                  numRecoBins,recoBinsArray);
00461     fRecoEnergyAllEvents_FD->Sumw2();
00462 
00463     //Checklist histos (CC numubars):
00464     fTrueEnergyCCOnlyEventsPQ_ND=new TH1D
00465                 ("TrueEnergyCCOnlyEventsPQ_ND",
00466                  "NuMuBar CC Selected Events with True Energy (NearDet)",
00467                  numTrueBins,trueBinsArray);
00468     fTrueEnergyCCOnlyEventsPQ_ND->Sumw2();
00469 
00470     fRecoEnergyCCOnlyEventsPQ_ND=new TH1D
00471                 ("RecoEnergyCCOnlyEventsPQ_ND",
00472                  "NuMuBar CC Selected Events with Reco Energy (NearDet)",
00473                  numRecoBins,recoBinsArray);
00474     fRecoEnergyCCOnlyEventsPQ_ND->Sumw2();
00475 
00476     fRecoEnergyAllEventsPQ_ND=new TH1D
00477                 ("RecoEnergyAllEventsPQ_ND",
00478                  "Selected Events with Reco Energy (NearDet)",
00479                  numRecoBins,recoBinsArray);
00480     fRecoEnergyAllEventsPQ_ND->Sumw2();
00481 
00482     fTrueEnergyCCOnlyEventsPQ_FD=new TH1D
00483                 ("TrueEnergyCCOnlyEventsPQ_FD",
00484                  "NuMuBar CC Selected Events with True Energy (FarDet)",
00485                  numTrueBins,trueBinsArray);
00486     fTrueEnergyCCOnlyEventsPQ_FD->Sumw2();
00487 
00488     fRecoEnergyCCOnlyEventsPQ_FD=new TH1D
00489                 ("RecoEnergyCCOnlyEventsPQ_FD",
00490                  "NuMuBar CC Selected Events with Reco Energy (FarDet)",
00491                  numRecoBins,recoBinsArray);
00492     fRecoEnergyCCOnlyEventsPQ_FD->Sumw2();
00493 
00494     fRecoEnergyAllEventsPQ_FD=new TH1D
00495                 ("RecoEnergyAllEventsPQ_FD",
00496                  "Selected Events with Reco Energy (FarDet)",
00497                  numRecoBins,recoBinsArray);
00498     fRecoEnergyAllEventsPQ_FD->Sumw2();
00499 
00500     //Checklist histos (CC numus + CC numubars):
00501     fTrueEnergyCCOnlyEventsAll_ND=new TH1D
00502                 ("TrueEnergyCCOnlyEventsAll_ND",
00503                  "All CC Selected Events with True Energy (NearDet)",
00504                  numTrueBins,trueBinsArray);
00505     fTrueEnergyCCOnlyEventsAll_ND->Sumw2();
00506 
00507     fRecoEnergyCCOnlyEventsAll_ND=new TH1D
00508                 ("RecoEnergyCCOnlyEventsAll_ND",
00509                  "All CC Selected Events with Reco Energy (NearDet)",
00510                  numRecoBins,recoBinsArray);
00511     fRecoEnergyCCOnlyEventsAll_ND->Sumw2();
00512 
00513     fRecoEnergyAllEventsAll_ND=new TH1D
00514                 ("RecoEnergyAllEventsAll_ND",
00515                  "Selected Events with Reco Energy (NearDet)",
00516                  numRecoBins,recoBinsArray);
00517     fRecoEnergyAllEventsAll_ND->Sumw2();
00518 
00519     fTrueEnergyCCOnlyEventsAll_FD=new TH1D
00520                 ("TrueEnergyCCOnlyEventsAll_FD",
00521                  "All CC Selected Events with True Energy (FarDet)",
00522                  numTrueBins,trueBinsArray);
00523     fTrueEnergyCCOnlyEventsAll_FD->Sumw2();
00524 
00525     fRecoEnergyCCOnlyEventsAll_FD=new TH1D
00526                 ("RecoEnergyCCOnlyEventsAll_FD",
00527                  "All CC Selected Events with Reco Energy (FarDet)",
00528                  numRecoBins,recoBinsArray);
00529     fRecoEnergyCCOnlyEventsAll_FD->Sumw2();
00530 
00531     fRecoEnergyAllEventsAll_FD=new TH1D
00532                 ("RecoEnergyAllEventsAll_FD",
00533                  "Selected Events with Reco Energy (FarDet)",
00534                  numRecoBins,recoBinsArray);
00535     fRecoEnergyAllEventsAll_FD->Sumw2();
00536 
00537     //Checklist histos (CC nutaus)
00538     fTrueEnergyCCOnlyEventsTau_FD=new TH1D
00539                 ("TrueEnergyCCOnlyEventsTau_FD",
00540                  "NuTau CC Selected Events with True Energy (FarDet)",
00541                  numTrueBins,trueBinsArray);
00542     fTrueEnergyCCOnlyEventsTau_FD->Sumw2();
00543 
00544     //Checklist histos (CC nutaubars)
00545     fTrueEnergyCCOnlyEventsTauPQ_FD=new TH1D
00546                 ("TrueEnergyCCOnlyEventsTauPQ_FD",
00547                  "NuTauBar CC Selected Events with True Energy (FarDet)",
00548                  numTrueBins,trueBinsArray);
00549     fTrueEnergyCCOnlyEventsTauPQ_FD->Sumw2();
00550 
00551     //Checklist histos (CC nutaus + nutaubars)
00552     fTrueEnergyCCOnlyEventsTauAll_FD=new TH1D
00553                 ("TrueEnergyCCOnlyEventsTauAll_FD",
00554                  "NuTau CC Selected Events with True Energy (FarDet)",
00555                  numTrueBins,trueBinsArray);
00556     fTrueEnergyCCOnlyEventsTauAll_FD->Sumw2();
00557 
00558     //DATA histos (CC numus):
00559     RecoEnergy_ND=new TH1D
00560                 ("RecoEnergy_ND","Reconstructed Energy - Data (NearDet)",
00561                  numRecoBins,recoBinsArray);
00562     RecoEnergy_FD=new TH1D
00563                 ("RecoEnergy_FD","Reconstructed Energy - Data (FarDet)",
00564                  numRecoBins,recoBinsArray);
00565 
00566     //DATA histos (CC numubars):
00567     RecoEnergyPQ_ND=new TH1D
00568                 ("RecoEnergyPQ_ND","Reconstructed Energy - Data (NearDet)",
00569                  numRecoBins,recoBinsArray);
00570     RecoEnergyPQ_FD=new TH1D
00571       ("RecoEnergyPQ_FD","Reconstructed Energy - Data (FarDet)",
00572        numRecoBins,recoBinsArray);
00573 
00574     //DATA histos (Pure Reco Energys)
00575     RecoEnergyTrueTauPQ_FD=new TH1D
00576       ("RecoEnergyTrueTauPQ_FD","Reconstructed Energy - Data (True anti tau FarDet)",
00577        numRecoBins,recoBinsArray);
00578     RecoEnergyTrueMuonPQ_FD=new TH1D
00579       ("RecoEnergyTrueMuonPQ_FD","Reconstructed Energy - Data (True anti muon FarDet)",
00580        numRecoBins,recoBinsArray);
00581     RecoEnergyTrueTauNQ_FD=new TH1D
00582   ("RecoEnergyTrueTauNQ_FD","Reconstructed Energy - Data (True tau FarDet)",
00583    numRecoBins,recoBinsArray);
00584     RecoEnergyTrueMuonNQ_FD=new TH1D
00585       ("RecoEnergyTrueMuonNQ_FD","Reconstructed Energy - Data (True muon FarDet)",
00586        numRecoBins,recoBinsArray);
00587 
00588     RecoEnergyTrueTauPQ_ND=new TH1D
00589   ("RecoEnergyTrueTauPQ_ND","Reconstructed Energy - Data (True anti tau NearDet)",
00590    numRecoBins,recoBinsArray);
00591     RecoEnergyTrueMuonPQ_ND=new TH1D
00592       ("RecoEnergyTrueMuonPQ_ND","Reconstructed Energy - Data (True anti muon NearDet)",
00593        numRecoBins,recoBinsArray);
00594     RecoEnergyTrueTauNQ_ND=new TH1D
00595       ("RecoEnergyTrueTauNQ_ND","Reconstructed Energy - Data (True tau NearDet)",
00596        numRecoBins,recoBinsArray);
00597     RecoEnergyTrueMuonNQ_ND=new TH1D
00598       ("RecoEnergyTrueMuonNQ_ND","Reconstructed Energy - Data (True muon NearDet)",
00599        numRecoBins,recoBinsArray);
00600 
00601     RecoEnergyWSTauPQ_FD=new TH1D
00602       ("RecoEnergyWSTauPQ_FD","Reconstructed Energy - Data (WS tau FarDet)",
00603        numRecoBins,recoBinsArray);
00604     RecoEnergyWSMuonPQ_FD=new TH1D
00605       ("RecoEnergyWSMuonPQ_FD","Reconstructed Energy - Data (WS muon FarDet)",
00606        numRecoBins,recoBinsArray);
00607     RecoEnergyWSTauNQ_FD=new TH1D
00608       ("RecoEnergyWSTauNQ_FD","Reconstructed Energy - Data (WS anti tau FarDet)",
00609        numRecoBins,recoBinsArray);
00610     RecoEnergyWSMuonNQ_FD=new TH1D
00611       ("RecoEnergyWSMuonNQ_FD","Reconstructed Energy - Data (WS anti muon FarDet)",
00612        numRecoBins,recoBinsArray);
00613 
00614     RecoEnergyWSTauPQ_ND=new TH1D
00615       ("RecoEnergyWSTauPQ_ND","Reconstructed Energy - Data (WS tau NearDet)",
00616        numRecoBins,recoBinsArray);
00617     RecoEnergyWSMuonPQ_ND=new TH1D
00618       ("RecoEnergyWSMuonPQ_ND","Reconstructed Energy - Data (WS muon NearDet)",
00619        numRecoBins,recoBinsArray);
00620     RecoEnergyWSTauNQ_FD=new TH1D
00621       ("RecoEnergyWSTauNQ_ND","Reconstructed Energy - Data (WS anti tau NearDet)",
00622        numRecoBins,recoBinsArray);
00623     RecoEnergyWSMuonNQ_ND=new TH1D
00624       ("RecoEnergyWSMuonNQ_ND","Reconstructed Energy - Data (WS anti muon NearDet)",
00625        numRecoBins,recoBinsArray);
00626 
00627 
00628     //DATA histos (CC numus + CC numubars):
00629     RecoEnergyAll_ND=new TH1D
00630                 ("RecoEnergyAll_ND","Reconstructed Energy - Data (NearDet)",
00631                  numRecoBins,recoBinsArray);
00632     RecoEnergyAll_FD=new TH1D
00633                 ("RecoEnergyAll_FD","Reconstructed Energy - Data (FarDet)",
00634                  numRecoBins,recoBinsArray);
00635 
00636         poKinPQ = new TH1F("poKinPQ", "poKin Value for PQ Events;poKIN", 100, -1., 1.);
00637         roIdPQ = new TH1F("roIdPQ", "roID Value for PQ Events;roID", 100, 0., 1.);
00638         dpIdPQ = new TH1F("dpIdPQ", "dpID Value for PQ Events;dpID", 100, -1., 1.);
00639         trkLPQ = new TH1F("trkLPQ", "trkLength for PQ Events;Track Length", 500, -10., 490.);
00640         relAPQ = new TH1F("relAPQ", "|Relative Angle - #pi| for PQ Events", 100, 0., 3.14159);
00641         relAPQ1 = new TH1F("relAPQ1","Relative Angle for PQ Events", 60, 0, 6.28);
00642         relAngVsTrkEndX = new TH2F("relAngVsTrkEndX", "Rel. Angle vs xTrkEnd for PQ", 100, -5., 5., 100, 2.0, 3.14159);
00643         relAngVsTrkEndY = new TH2F("relAngVsTrkEndY", "Rel. Angle vs yTrkEnd for PQ", 100, -3.5, 3.5, 100, 2.0, 3.14159);
00644         relAngVsTrkEndZ = new TH2F("relAngVsTrkEndZ", "Rel. Angle vs zTrkEnd for PQ", 100, 0., 18., 100, 2.0, 3.14159);
00645         relANQ = new TH1F("relANQ", "|Relative Angle - #pi| for NQ Events", 100, 0., 3.14159);
00646         qp_sigmaqpPQ = new TH1F("qp_sigmaqpPQ", "qp/#sigma qp for PQ Events", 320, -16, 16.);
00647         majCPQ = new TH1F("majCPQ", "Majority Curvature for PQ Events;Majority Curvature", 100, -1., 1.);
00648         probPQ = new TH1F("probPQ", "prob for PQ Events;Prob", 100, 0., 1.);
00649   }
00650 
00651   Float_t weight=nu.rw;
00652 
00653   if (nu.simFlag==SimFlag::kMC) {
00654     if (nu.charge==-1) {
00655       //fill histos
00656       if (nu.detector==Detector::kNear) {
00657         fRecoEnergyAllEvents_ND->Fill(nu.energy,weight);
00658         if(nu.iaction==1 && nu.inu==14) {
00659           fRecoEnergyCCOnlyEvents_ND->Fill(nu.energy,weight);
00660           fTrueEnergyCCOnlyEvents_ND->Fill(fabs(nu.energyMC),weight);
00661           fPurity_ND->Fill(nu.energy,weight);
00662           fEfficiency_ND->Fill(fabs(nu.energyMC),weight);
00663           fRecoVsTrueEnergy_ND->Fill(fabs(nu.energyMC),
00664                                      nu.energy,weight);
00665         }
00666         else if (nu.iaction==1 && nu.inu==-14){
00667           fCCContamination_ND->Fill(nu.energy,weight);
00668         }
00669         else {
00670           fNCContamination_ND->Fill(nu.energy,weight);
00671         }
00672       }
00673       else if (nu.detector==Detector::kFar) {
00674         fRecoEnergyAllEvents_FD->Fill(nu.energy,weight);
00675         if(nu.iaction==1 && nu.inu==14) {
00676           fRecoEnergyCCOnlyEvents_FD->Fill(nu.energy,weight);
00677           fTrueEnergyCCOnlyEvents_FD->Fill(fabs(nu.energyMC),weight);
00678           fPurity_FD->Fill(nu.energy,weight);
00679           fEfficiency_FD->Fill(fabs(nu.energyMC),weight);
00680           fRecoVsTrueEnergy_FD->Fill(fabs(nu.energyMC),
00681                                      nu.energy,weight);
00682         }
00683         else if (nu.iaction==1 && nu.inu==-14){
00684           fCCContamination_FD->Fill(fabs(nu.energyMC),weight);
00685           fCCContaminationRecoVsTrue_FD->Fill(fabs(nu.energyMC),
00686                                               nu.energy,
00687                                               weight);
00688         }
00689         else if (nu.iaction==1 && nu.inu==16 && nu.inunoosc==14){
00690           fEfficiencyTau_FD->Fill(fabs(nu.energyMC),weight);
00691           fTrueEnergyCCOnlyEventsTau_FD->Fill(fabs(nu.energyMC),weight);
00692           fRecoVsTrueEnergyTau_FD->Fill(fabs(nu.energyMC),
00693                                         nu.energy,weight);
00694         }
00695         else {
00696           fNCContamination_FD->Fill(nu.energy,weight);
00697         }
00698       }
00699       else cout<<"ahhhh, det="<<nu.detector<<endl;
00700     }
00701     if (nu.charge==1) {
00702       //fill histos
00703       if (nu.detector==Detector::kNear) {
00704         fRecoEnergyAllEventsPQ_ND->Fill(nu.energy,weight);
00705         if(nu.iaction==1 && nu.inu==-14) {
00706           fRecoEnergyCCOnlyEventsPQ_ND->Fill(nu.energy,weight);
00707           fTrueEnergyCCOnlyEventsPQ_ND->Fill(fabs(nu.energyMC),weight);
00708           fPurityPQ_ND->Fill(nu.energy,weight);
00709           fEfficiencyPQ_ND->Fill(fabs(nu.energyMC),weight);
00710           fRecoVsTrueEnergyPQ_ND->Fill(fabs(nu.energyMC),
00711                                        nu.energy,weight);
00712         }
00713         else if (nu.iaction==1 && nu.inu==14){
00714           fCCContaminationPQ_ND->Fill(nu.energy,weight);
00715         }
00716         else {
00717           fNCContaminationPQ_ND->Fill(nu.energy,weight);
00718         }
00719       }
00720       else if (nu.detector==Detector::kFar) {
00721         fRecoEnergyAllEventsPQ_FD->Fill(nu.energy,weight);
00722         if(nu.iaction==1 && nu.inu==-14) {
00723           fRecoEnergyCCOnlyEventsPQ_FD->Fill(nu.energy,weight);
00724           fTrueEnergyCCOnlyEventsPQ_FD->Fill(fabs(nu.energyMC),weight);
00725           fPurityPQ_FD->Fill(nu.energy,weight);
00726           fEfficiencyPQ_FD->Fill(fabs(nu.energyMC),weight);
00727           fRecoVsTrueEnergyPQ_FD->Fill(fabs(nu.energyMC),
00728                                        nu.energy,weight);
00729         }
00730         else if (nu.iaction==1 && nu.inu==14){
00731           fCCContaminationPQ_FD->Fill(fabs(nu.energyMC),weight);
00732           fCCContaminationRecoVsTruePQ_FD->Fill(fabs(nu.energyMC),
00733                                                 nu.energy,
00734                                                 weight);
00735         }
00736         else if (nu.iaction==1 && nu.inu==-16 && nu.inunoosc==-14){
00737           fEfficiencyTauPQ_FD->Fill(fabs(nu.energyMC),weight);
00738           fTrueEnergyCCOnlyEventsTauPQ_FD->Fill(fabs(nu.energyMC),weight);
00739           fRecoVsTrueEnergyTauPQ_FD->Fill(fabs(nu.energyMC),
00740                                           nu.energy,weight);
00741         }
00742         else {
00743           fNCContaminationPQ_FD->Fill(nu.energy,weight);
00744         }
00745       }
00746       else cout<<"ahhhh, det="<<nu.detector<<endl;
00747     }
00748     if ((nu.charge==1) || (-1==nu.charge)) {
00749       //fill histos
00750       if (nu.detector==Detector::kNear) {
00751         fRecoEnergyAllEventsAll_ND->Fill(nu.energy,weight);
00752         if(nu.iaction==1 && (nu.inu==-14 || nu.inu == 14)) {
00753           fRecoEnergyCCOnlyEventsAll_ND->Fill(nu.energy,weight);
00754           fTrueEnergyCCOnlyEventsAll_ND->Fill(fabs(nu.energyMC),weight);
00755           fPurityAll_ND->Fill(nu.energy,weight);
00756           fEfficiencyAll_ND->Fill(fabs(nu.energyMC),weight);
00757           fRecoVsTrueEnergyAll_ND->Fill(fabs(nu.energyMC),
00758                                         nu.energy,weight);
00759         }
00760       }
00761       else if (nu.detector==Detector::kFar) {
00762         fRecoEnergyAllEventsAll_FD->Fill(nu.energy,weight);
00763         if(nu.iaction==1 && (nu.inu==-14 || nu.inu == 14)) {
00764           fRecoEnergyCCOnlyEventsAll_FD->Fill(nu.energy,weight);
00765           fTrueEnergyCCOnlyEventsAll_FD->Fill(fabs(nu.energyMC),weight);
00766           fPurityAll_FD->Fill(nu.energy,weight);
00767           fEfficiencyAll_FD->Fill(fabs(nu.energyMC),weight);
00768           fRecoVsTrueEnergyAll_FD->Fill(fabs(nu.energyMC),
00769                                         nu.energy,weight);
00770         }
00771         else if (nu.iaction==1 && (nu.inu==-16 || nu.inu == 16)
00772                  && (nu.inunoosc == -14 || nu.inunoosc == 14)){
00773           fEfficiencyTauAll_FD->Fill(fabs(nu.energyMC),weight);
00774           fTrueEnergyCCOnlyEventsTauAll_FD->Fill(fabs(nu.energyMC),weight);
00775           fRecoVsTrueEnergyTauAll_FD->Fill(fabs(nu.energyMC),
00776                                           nu.energy,weight);
00777         }
00778       }
00779       else cout<<"ahhhh, det="<<nu.detector<<endl;
00780     }
00781   }
00782 
00783   Bool_t oscillateData = false;
00784   Double_t oscWeight = 1.0;
00785   Double_t s2t = 1.0;
00786   Double_t dm2 = 2.3e-3;
00787   if (oscillateData && nu.energyMC && nu.iaction==1 && abs(nu.inu)==14){
00788     oscWeight = NuUtilities::OscillationWeight(nu.energyMC, dm2, s2t);
00789   }
00790 
00791   if (nu.charge==-1) {
00792     if (nu.detector==Detector::kNear) {
00793       RecoEnergy_ND->Fill(nu.energy,weight);
00794     }
00795     else if (nu.detector==Detector::kFar) {
00796       RecoEnergy_FD->Fill(nu.energy,weight*oscWeight);
00797     }
00798     else cout<<"ahhhh, det="<<nu.detector<<endl;
00799   }
00800   if (nu.charge==1) {
00801     if (nu.detector==Detector::kNear) {
00802       RecoEnergyPQ_ND->Fill(nu.energy,weight);
00803     }
00804     else if (nu.detector==Detector::kFar) {
00805       RecoEnergyPQ_FD->Fill(nu.energy,weight*oscWeight);
00806     }
00807     else cout<<"ahhhh, det="<<nu.detector<<endl;
00808   }
00809   if(nu.detector==Detector::kFar&&nu.iaction==1){
00810     if(nu.inu==-16&&nu.charge==1&& nu.inunoosc==-14){
00811       RecoEnergyTrueTauPQ_FD->Fill(nu.energy,weight*oscWeight);
00812     }
00813     if(nu.inu==16&&nu.charge==1&& nu.inunoosc==14){
00814       RecoEnergyWSTauPQ_FD->Fill(nu.energy,weight*oscWeight);
00815     }
00816     if(nu.inu==-14&&nu.charge==1){
00817       RecoEnergyTrueMuonPQ_FD->Fill(nu.energy,weight*oscWeight);
00818     }
00819     if(nu.inu==14&&nu.charge==1){
00820       RecoEnergyWSMuonPQ_FD->Fill(nu.energy,weight*oscWeight);
00821     }
00822     if(nu.inu==16&&nu.charge==-1&& nu.inunoosc==14){
00823       RecoEnergyTrueTauNQ_FD->Fill(nu.energy,weight*oscWeight);
00824     }
00825     if(nu.inu==-16&&nu.charge==-1&& nu.inunoosc==-14){
00826       RecoEnergyWSTauNQ_FD->Fill(nu.energy,weight*oscWeight);
00827     }
00828     if(nu.inu==14&&nu.charge==-1){
00829       RecoEnergyTrueMuonNQ_FD->Fill(nu.energy,weight*oscWeight);
00830     }
00831     if(nu.inu==-14&&nu.charge==-1){
00832       RecoEnergyWSMuonNQ_FD->Fill(nu.energy,weight*oscWeight);
00833     }
00834   }
00835   else if(nu.detector==Detector::kNear&&nu.iaction==1){
00836     if(nu.inu==-16&&nu.charge==1&& nu.inunoosc==-14){
00837       RecoEnergyTrueTauPQ_ND->Fill(nu.energy,weight*oscWeight);
00838     }
00839     if(nu.inu==16&&nu.charge==1&& nu.inunoosc==14){
00840       RecoEnergyWSTauPQ_ND->Fill(nu.energy,weight*oscWeight);
00841     }
00842     if(nu.inu==-14&&nu.charge==1){
00843       RecoEnergyTrueMuonPQ_ND->Fill(nu.energy,weight*oscWeight);
00844     }
00845     if(nu.inu==14&&nu.charge==1){
00846       RecoEnergyWSMuonPQ_ND->Fill(nu.energy,weight*oscWeight);
00847     }
00848     if(nu.inu==16&&nu.charge==-1&& nu.inunoosc==14){
00849       RecoEnergyTrueTauNQ_ND->Fill(nu.energy,weight*oscWeight);
00850     }
00851     if(nu.inu==-16&&nu.charge==-1&& nu.inunoosc==-14){
00852       RecoEnergyWSTauNQ_ND->Fill(nu.energy,weight*oscWeight);
00853     }
00854     if(nu.inu==14&&nu.charge==-1){
00855       RecoEnergyTrueMuonNQ_ND->Fill(nu.energy,weight*oscWeight);
00856     }
00857     if(nu.inu==-14&&nu.charge==-1){
00858       RecoEnergyWSMuonNQ_ND->Fill(nu.energy,weight*oscWeight);
00859     }
00860   }
00861 
00862   if (nu.charge==1 || -1==nu.charge) {
00863     if (nu.detector==Detector::kNear) {
00864       RecoEnergyAll_ND->Fill(nu.energy,weight);
00865     }
00866     else if (nu.detector==Detector::kFar) {
00867       RecoEnergyAll_FD->Fill(nu.energy,weight*oscWeight);
00868     }
00869     else cout<<"ahhhh, det="<<nu.detector<<endl;
00870   }
00871 
00872 
00873     // Do the cut variables...
00874     if (nu.charge == 1) {
00875       if (weight != 0.)
00876         {
00877           //      cout<<"check if weights are correct = "<<weight<<endl;
00878         }
00879     poKinPQ->Fill(nu.poIDKin,weight);
00880     roIdPQ->Fill(nu.roID,weight);
00881     dpIdPQ->Fill(nu.dpID,weight);
00882     trkLPQ->Fill(nu.trkLength,weight);
00883     relAPQ->Fill(abs(nu.relativeAngle-TMath::Pi()),weight);
00884     relAPQ1->Fill(nu.relativeAngle,weight);
00885     relAngVsTrkEndX->Fill(nu.xTrkEnd,abs(nu.relativeAngle-TMath::Pi()),weight);
00886     relAngVsTrkEndY->Fill(nu.yTrkEnd,abs(nu.relativeAngle-TMath::Pi()),weight);
00887     relAngVsTrkEndZ->Fill(nu.zTrkEnd,abs(nu.relativeAngle-TMath::Pi()),weight);
00888     qp_sigmaqpPQ->Fill(1./nu.sigqp_qp,weight);
00889     majCPQ->Fill(nu.smoothMajC,weight);
00890     probPQ->Fill(nu.prob,weight);
00891   }
00892 }
00893 
00894 //......................................................................
00895 
00896 void NuHistos::
00897 FillMatrixMethodNCHistos(const NuEvent& nu,
00898                          const NuBinningScheme::NuBinningScheme_t
00899                          binningScheme) const
00900 {
00901   MAXMSG("NuHistos",Msg::kInfo,1)
00902     <<"Running FillMatrixMethodNCHistos..."<<endl;
00903 
00904   //DATA histos (NCs)
00905   static TH1D* recoEnergy_ND_NC = 0;
00906   static TH1D* recoEnergy_FD_NC = 0;
00907   static TH2D* recoVsTrueEnergy_ND_NC = 0;
00908   static TH2D* recoVsTrueEnergy_FD_NC = 0;
00909 
00910   static TH2D* recoVsTrueEnergy_truly[kNumTruths][2] = {{0,},}; // ND, FD
00911 
00912   if(!recoEnergy_ND_NC){ // If not initialized
00913     MSG("NuHistos",Msg::kInfo)
00914       <<"Creating histograms on first run of FillMatrixMethodNCHistos..."
00915       <<" gDirectory gives:"<<endl;
00916           gDirectory->Print();
00917 
00918     const NuUtilities cuts;
00919     const vector<Double_t> recoBins = cuts.RecoBins(binningScheme);
00920     const Int_t numRecoBins = recoBins.size()-1;
00921     const Double_t* recoBinsArray = &recoBins[0];
00922 
00923     const vector<Double_t> trueBins = cuts.TrueBins(binningScheme);
00924     const Int_t numTrueBins = trueBins.size()-1;
00925     const Double_t* trueBinsArray = &trueBins[0];
00926 
00927     recoEnergy_ND_NC=new TH1D
00928       ("RecoEnergy_ND_NC","Reconstructed NC Energy - Data (NearDet)",
00929        numRecoBins, recoBinsArray);
00930     recoEnergy_FD_NC=new TH1D
00931       ("RecoEnergy_FD_NC","Reconstructed NC Energy - Data (FarDet)",
00932        numRecoBins, recoBinsArray);
00933 
00934     recoVsTrueEnergy_ND_NC=new TH2D
00935       ("RecoVsTrueEnergy_ND_NC",
00936        "Reco Vs True Energy (NearDet);True neutrino energy;Reconstructed energy",
00937        numTrueBins, trueBinsArray,
00938        numRecoBins, recoBinsArray);
00939     recoVsTrueEnergy_ND_NC->Sumw2();
00940 
00941     recoVsTrueEnergy_FD_NC=new TH2D
00942       ("RecoVsTrueEnergy_FD_NC",
00943        "Reco Vs True Energy (FarDet);True neutrino energy;Reconstructed energy",
00944        numTrueBins, trueBinsArray,
00945        numRecoBins, recoBinsArray);
00946     recoVsTrueEnergy_FD_NC->Sumw2();
00947 
00948     for(int truth = 0; truth < kNumTruths; ++truth){
00949       recoVsTrueEnergy_truly[truth][0] =
00950         new TH2D("RecoVsTrueEnergy_truly"+truthNames[truth]+"_ND_NC",
00951                  "Reco Vs True Energy (NearDet);"
00952                  "True neutrino energy;Reconstructed energy",
00953                  numTrueBins, trueBinsArray, numRecoBins, recoBinsArray);
00954 
00955       recoVsTrueEnergy_truly[truth][1] =
00956         new TH2D("RecoVsTrueEnergy_truly"+truthNames[truth]+"_FD_NC",
00957                  "Reco Vs True Energy (FarDet);"
00958                  "True neutrino energy; Reconstructed energy",
00959                  numTrueBins, trueBinsArray, numRecoBins, recoBinsArray);
00960     } // end for truth
00961   } // end if not initialized
00962 
00963   assert(nu.detector == Detector::kNear || nu.detector == Detector::kFar);
00964 
00965   // TODO CJB Use NC energy
00966   //  if(nu.detector == Detector::kNear) recoEnergy_ND_NC->Fill(nu.energy, nu.rw);
00967   //  if(nu.detector == Detector::kFar)  recoEnergy_FD_NC->Fill(nu.energy, nu.rw);
00968 
00969   ENCTruth_t truth = kNumTruths;
00970   if(abs(nu.inu) == 12 && nu.iaction == 1) truth = kCCe;
00971   if(nu.inu == +14 && nu.iaction == 1) truth = kCCmu;
00972   if(nu.inu == -14 && nu.iaction == 1) truth = kCCmubar;
00973   if(nu.inu == +16 && nu.iaction == 1) truth = kCCtau;
00974   if(nu.inu == -16 && nu.iaction == 1) truth = kCCtaubar;
00975   if(nu.inu > 0    && nu.iaction == 0) truth = kNC;
00976   if(nu.inu < 0    && nu.iaction == 0) truth = kNCbar;
00977   assert(truth != kNumTruths);
00978 
00979   const int det = int(nu.detector == Detector::kFar);
00980 
00981   // TODO CJB Use NC energy
00982   recoVsTrueEnergy_truly[truth][det]->Fill(fabs(nu.neuEnMC), nu.energy, nu.rw);
00983 
00984 
00985   switch(nu.detector){
00986   case Detector::kNear:
00987     // TODO CJB use NC energy
00988     //    recoEnergy_ND_NC->Fill(nu.energyNC, nu.rw);
00989     //    recoVsTrueEnergy_ND_NC->Fill(fabs(nu.neuEnMC), nu.energyNC, nu.rw);
00990     recoEnergy_ND_NC->Fill(nu.energy, nu.rw);
00991     recoVsTrueEnergy_ND_NC->Fill(fabs(nu.neuEnMC), nu.energy, nu.rw);
00992     break;
00993   case Detector::kFar:
00994     // TODO CJB use NC energy
00995     //    recoEnergy_FD_NC->Fill(nu.energyNC, nu.rw);
00996     //    recoVsTrueEnergy_FD_NC->Fill(fabs(nu.neuEnMC), nu.energyNC, nu.rw);
00997     recoEnergy_FD_NC->Fill(nu.energy, nu.rw);
00998     recoVsTrueEnergy_FD_NC->Fill(fabs(nu.neuEnMC), nu.energy, nu.rw);
00999     break;
01000   default:
01001     cout<<"ahhhh, det="<<nu.detector<<endl;
01002     assert(0 && "bad detector");
01003   }
01004 }
01005 
01006 //......................................................................
01007 
01008 
01009 
01010 
01011 /*
01012 
01013   //Helper histos:
01014   TH2D *fRecoVsTrueEnergy_ND; // conversion from reco->true
01015   TH2D *fRecoVsTrueEnergy_FD; // for ND/FD
01016   TH2D *fFDVsNDMatrix;        // ND->FD number of events conversion
01017   TH2D *fFDVsNDMatrixRW;      // ND->FD number of events conversion+new nearwt
01018   TH2D *fFDVsNDMatrixXSec;    // ND->FD number of events conversion with xsec folded in
01019   TH2D *fFDVsNDMatrixXSecRW;  // ND->FD number of events conversion with xsec + new nearwt
01020   TH1D *fEfficiency_ND;       // FD and ND efficiency (with true energy)
01021   TH1D *fEfficiency_FD;
01022   TH1D *fPurity_ND;           // FD and ND purity (with reco energy)
01023   TH1D *fPurity_FD;
01024   TH1D *fXSec_CC;             // numu CC cross-section with energy
01025   TH1D *fFracErrOnPred;       // fractional error on predicted FD spectrum with energy
01026 
01027   //Check list histos:
01028   TH1D *fRecoEnergyAllEvents_ND;
01029   TH1D *fRecoEnergyCCOnlyEvents_ND;
01030   TH1D *fTrueEnergyCCOnlyEvents_ND;
01031   TH1D *fTrueEnergyTrueCCFidEvents_ND;
01032   TH1D *fTrueEnergyNuFlux_ND;
01033   TH1D *fTrueEnergyNuFluxRW_ND;
01034   TH1D *fTrueEnergyNuFlux_FD;
01035   TH1D *fTrueEnergyNuFluxRW_FD;
01036   TH1D *fTrueEnergyCCFlux_ND;
01037   TH1D *fTrueEnergyCCFluxRW_ND;
01038   TH1D *fTrueEnergyCCFlux_FD;
01039   TH1D *fTrueEnergyCCFluxRW_FD;
01040   TH1D *fTrueEnergyTrueCCFidEvents_FD;
01041   TH1D *fTrueEnergyCCOnlyEvents_FD;
01042   TH1D *fRecoEnergyCCOnlyEvents_FD;
01043   TH1D *fRecoEnergyAllEvents_FD;
01044 
01045 
01046     //helper histos:
01047     fRecoVsTrueEnergy_ND=new TH2D
01048       ("RecoVsTrueEnergy_ND",
01049        "Reco Vs True Energy (NearDet)",
01050        fNBins,fXBins,fNBins,fXBins);
01051     fRecoVsTrueEnergy_ND->Sumw2();
01052 
01053     fRecoVsTrueEnergy_FD=new TH2D
01054       ("RecoVsTrueEnergy_FD",
01055        "Reco Vs True Energy (FarDet)",
01056        fNBins,fXBins,fNBins,fXBins);
01057     fRecoVsTrueEnergy_FD->Sumw2();
01058 
01059     fEfficiency_ND=new TH1D
01060       ("Efficiency_ND",
01061        "NuMu CC Selection Efficiency with True Energy (NearDet)",
01062        fNBins,fXBins);
01063     fEfficiency_ND->Sumw2();
01064 
01065     fEfficiency_FD=new TH1D
01066       ("Efficiency_FD",
01067        "NuMu CC Selection Efficiency with True Energy (FarDet)",
01068        fNBins,fXBins);
01069     fEfficiency_FD->Sumw2();
01070 
01071     fPurity_ND=new TH1D
01072       ("Purity_ND",
01073        "NuMu CC Selection Purity with True Energy (NearDet)",
01074        fNBins,fXBins);
01075     fPurity_ND->Sumw2();
01076 
01077     fPurity_FD=new TH1D
01078       ("Purity_FD",
01079        "NuMu CC Selection Purity with True Energy (FarDet)",
01080        fNBins,fXBins);
01081     fPurity_FD->Sumw2();
01082 
01083     fFDVsNDMatrix=new TH2D
01084       ("FDVsNDMatrix",
01085        "Number of FD Vs ND Events with True Energy",
01086        fNBins,fXBins,fNBins,fXBins);
01087     fFDVsNDMatrix->Sumw2();
01088 
01089     fFDVsNDMatrixRW=new TH2D
01090       ("FDVsNDMatrixRW",
01091        "Number of FD Vs ND Events with True Energy (with Near Reweight)",
01092        fNBins,fXBins,fNBins,fXBins);
01093     fFDVsNDMatrixRW->Sumw2();
01094 
01095     fFDVsNDMatrixXSec=new TH2D
01096       ("FDVsNDMatrixXSec",
01097        "Number of FD Vs ND Events with True Energy (with XSec)",
01098        fNBins,fXBins,fNBins,fXBins);
01099     fFDVsNDMatrixXSec->Sumw2();
01100 
01101     fFDVsNDMatrixXSecRW=new TH2D
01102       ("FDVsNDMatrixXSecRW",
01103        "Number of FD Vs ND Events with True Energy (with XSec + Near Reweight)",
01104        fNBins,fXBins,fNBins,fXBins);
01105     fFDVsNDMatrixXSecRW->Sumw2();
01106 
01107     fXSec_CC=new TH1D
01108       ("XSec_CC","NuMu CC XSection with True Energy",
01109        fNBins,fXBins);
01110     fXSec_CC->Sumw2();
01111 
01112     fFracErrOnPred=new TH1D
01113       ("FracErrOnPred",
01114        "Fractional Error on Energy Spectrum with Reco Energy",
01115        fNBins,fXBins);
01116     fFracErrOnPred->Sumw2();
01117 
01118     //Check list histos:
01119     fRecoEnergyAllEvents_ND=new TH1D
01120       ("RecoEnergyAllEvents_ND",
01121        "Selected Events with Reco Energy (NearDet)",
01122        fNBins,fXBins);
01123     fRecoEnergyAllEvents_ND->Sumw2();
01124 
01125     fRecoEnergyCCOnlyEvents_ND=new TH1D
01126       ("RecoEnergyCCOnlyEvents_ND",
01127        "NuMu CC Selected Events with Reco Energy (NearDet)",
01128        fNBins,fXBins);
01129     fRecoEnergyCCOnlyEvents_ND->Sumw2();
01130 
01131     fTrueEnergyCCOnlyEvents_ND=new TH1D
01132       ("TrueEnergyCCOnlyEvents_ND",
01133        "NuMu CC Selected Events with True Energy (NearDet)",
01134        fNBins,fXBins);
01135     fTrueEnergyCCOnlyEvents_ND->Sumw2();
01136 
01137     fTrueEnergyTrueCCFidEvents_ND=new TH1D
01138       ("TrueEnergyTrueCCFidEvents_ND",
01139        "True Fid NuMu CC Events with True Energy (NearDet)",
01140        fNBins,fXBins);
01141     fTrueEnergyTrueCCFidEvents_ND->Sumw2();
01142 
01143     fTrueEnergyNuFlux_ND=new TH1D
01144       ("TrueEnergyNuFlux_ND",
01145        "Neutrino Flux with True Energy (NearDet)",
01146        fNBins,fXBins);
01147     fTrueEnergyNuFlux_ND->Sumw2();
01148 
01149     fTrueEnergyNuFluxRW_ND=new TH1D
01150       ("TrueEnergyNuFluxRW_ND",
01151        "Neutrino Flux with True Energy (NearDet with Reweighting)",
01152        fNBins,fXBins);
01153     fTrueEnergyNuFluxRW_ND->Sumw2();
01154 
01155     fTrueEnergyNuFlux_FD=new TH1D
01156       ("TrueEnergyNuFlux_FD",
01157        "Neutrino Flux with True Energy (FarDet)",
01158        fNBins,fXBins);
01159     fTrueEnergyNuFlux_FD->Sumw2();
01160 
01161     fTrueEnergyNuFluxRW_FD=new TH1D
01162       ("TrueEnergyNuFluxRW_FD",
01163        "Neutrino Flux with True Energy (FarDet with Reweighting)",
01164        fNBins,fXBins);
01165     fTrueEnergyNuFluxRW_FD->Sumw2();
01166 
01167     fTrueEnergyCCFlux_ND=new TH1D
01168       ("TrueEnergyCCFlux_ND",
01169        "NuMu CC Flux with True Energy (NearDet)",
01170        fNBins,fXBins);
01171     fTrueEnergyCCFlux_ND->Sumw2();
01172 
01173     fTrueEnergyCCFluxRW_ND=new TH1D
01174       ("TrueEnergyCCFluxRW_ND",
01175        "NuMu CC Flux with True Energy with (NearDet with Reweighting)",
01176        fNBins,fXBins);
01177     fTrueEnergyCCFluxRW_ND->Sumw2();
01178 
01179     fTrueEnergyCCFlux_FD=new TH1D
01180       ("TrueEnergyCCFlux_FD",
01181        "NuMu CC Flux with True Energy (FarDet)",
01182        fNBins,fXBins);
01183     fTrueEnergyCCFlux_FD->Sumw2();
01184 
01185     fTrueEnergyCCFluxRW_FD=new TH1D
01186       ("TrueEnergyCCFluxRW_FD",
01187        "NuMu CC Flux with True Energy (FarDet with Reweighing)",
01188        fNBins,fXBins);
01189     fTrueEnergyCCFluxRW_FD->Sumw2();
01190 
01191     fTrueEnergyTrueCCFidEvents_FD=new TH1D
01192       ("TrueEnergyTrueCCFidEvents_FD",
01193        "True Fid NuMu CC Events with True Energy (FarDet)",
01194        fNBins,fXBins);
01195     fTrueEnergyTrueCCFidEvents_FD->Sumw2();
01196       fTrueEnergyCCOnlyEvents_FD=new TH1D
01197       ("TrueEnergyCCOnlyEvents_FD",
01198        "NuMu CC Selected Events with True Energy (FarDet)",
01199        fNBins,fXBins);
01200     fTrueEnergyCCOnlyEvents_FD->Sumw2();
01201 
01202     fRecoEnergyCCOnlyEvents_FD=new TH1D
01203       ("RecoEnergyCCOnlyEvents_FD",
01204        "NuMu CC Selected Events with Reco Energy (FarDet)",
01205        fNBins,fXBins);
01206     fRecoEnergyCCOnlyEvents_FD->Sumw2();
01207 
01208     fRecoEnergyAllEvents_FD=new TH1D
01209       ("RecoEnergyAllEvents_FD",
01210        "Selected Events with Reco Energy (FarDet)",
01211        fNBins,fXBins);
01212     fRecoEnergyAllEvents_FD->Sumw2();
01213 
01214 
01215 
01216   static TH1D* fTrueEnergyNuFlux_FD=0;
01217   static TH1D* fTrueEnergyNuFluxRW_FD=0;
01218   static TH1D* fTrueEnergyCCFlux_FD=0;
01219   static TH1D* fTrueEnergyCCFluxRW_FD=0;
01220   static TH1D* fTrueEnergyTrueCCFidEvents_FD=0;
01221 
01222 
01223   //Helper histos:
01224   static TH2D* fRecoVsTrueEnergy_ND=0;//conversion from reco->true
01225   static TH1D* fEfficiency_ND=0;//efficiency (with true energy)
01226   static TH1D* fPurity_ND=0;//purity (with reco energy)
01227   static TH1D* fXSec_CC=0;//numu CC cross-section with energy
01228   static TH1D* fFracErrOnPred=0;//fract. error on predicted FD spectrum with energy
01229   static TH2D* fRecoVsTrueEnergy_ND=0; // conversion from reco->true
01230   static TH2D* fRecoVsTrueEnergy_FD=0; // for ND/FD
01231   static TH2D* fFDVsNDMatrix=0;        // ND->FD number of events conversion
01232   static TH2D* fFDVsNDMatrixRW=0;      // ND->FD number of events conversion+new nearwt
01233   static TH2D* fFDVsNDMatrixXSec=0;    // ND->FD number of events conversion with xsec folded in
01234   static TH2D* fFDVsNDMatrixXSecRW=0;  // ND->FD number of events conversion with xsec + new nearwt
01235   static TH1D* fEfficiency_ND=0;       // FD and ND efficiency (with true energy)
01236   static TH1D* fEfficiency_FD=0;
01237   static TH1D* fPurity_ND=0;           // FD and ND purity (with reco energy)
01238   static TH1D* fPurity_FD=0;
01239   static TH1D* fXSec_CC=0;             // numu CC cross-section with energy
01240   static TH1D* fFracErrOnPred=0;       // fractional error on predicted FD spectrum with energy
01241 
01242   //Check list histos:
01243   static TH1D* fRecoEnergyAllEvents_ND=0;
01244   static TH1D* fRecoEnergyCCOnlyEvents_ND=0;
01245   static TH1D* fTrueEnergyCCOnlyEvents_ND=0;
01246   static TH1D* fTrueEnergyTrueCCFidEvents_ND=0;
01247   static TH1D* fTrueEnergyNuFlux_ND=0;
01248   static TH1D* fTrueEnergyNuFluxRW_ND=0;
01249   static TH1D* fTrueEnergyNuFlux_FD=0;
01250   static TH1D* fTrueEnergyNuFluxRW_FD=0;
01251   static TH1D* fTrueEnergyCCFlux_ND=0;
01252   static TH1D* fTrueEnergyCCFluxRW_ND=0;
01253   static TH1D* fTrueEnergyCCFlux_FD=0;
01254   static TH1D* fTrueEnergyCCFluxRW_FD=0;
01255   static TH1D* fTrueEnergyTrueCCFidEvents_FD=0;
01256   static TH1D* fTrueEnergyCCOnlyEvents_FD=0;
01257   static TH1D* fRecoEnergyCCOnlyEvents_FD=0;
01258   static TH1D* fRecoEnergyAllEvents_FD=0;
01259 
01260 
01261 */
01262 
01263 //......................................................................
01264 
01265 Float_t NuHistos::CalcPOTsFromHistos(std::string sNamePostfix) const
01266 {
01267   //Method to Calc POTs from histograms
01268 
01269   //string to use for histo names
01270   string s="";
01271 
01272 
01273   //get the detector
01274   s="hDetector";
01275   s+=sNamePostfix;
01276   TH1F* hDetector=(TH1F*)gROOT->FindObjectAny(s.c_str());
01277   Int_t det=-1;
01278   if (hDetector) det=static_cast<Int_t>(hDetector->GetMean());
01279   else cout<<"Ahhh, no hDetector"<<endl;
01280 
01281   //get the simflag
01282   s="hSimFlag";
01283   s+=sNamePostfix;
01284   TH1F* hSimFlag=(TH1F*)gROOT->FindObjectAny(s.c_str());
01285   Int_t simFlag=-1;
01286   if (hSimFlag) simFlag=static_cast<Int_t>(hSimFlag->GetMean());
01287   else cout<<"Ahhh, no hSimFlag"<<endl;
01288 
01289   cout<<"POT Counting: found Detector="<<det
01290       <<", SimFlag="<<simFlag<<endl;
01291 
01292   //get the pots
01293   s="hPottrtgtd";
01294   s+=sNamePostfix;
01295   TH1F* hPot=(TH1F*)gROOT->FindObjectAny(s.c_str());
01296   Float_t avPot=-1;
01297   if (hPot) avPot=hPot->GetMean();
01298   else cout<<"Ahhh, no hPottrtgtd histogram"<<endl;
01299 
01300   //get the pots that failed the cuts
01301   s="hPotBadtrtgtd";
01302   s+=sNamePostfix;
01303   TH1F* hPotBad=(TH1F*)gROOT->FindObjectAny(s.c_str());
01304   Float_t avPotBad=-1;
01305   if (hPotBad) avPotBad=hPotBad->GetMean();
01306   else cout<<"Ahhh, no hPotBadtrtgtd histogram"<<endl;
01307     
01308   //get the number of runs/subruns
01309   s="hRun";
01310   s+=sNamePostfix;
01311   TH1F* hRun=(TH1F*)gROOT->FindObjectAny(s.c_str());
01312   Float_t runs=-1;
01313   if (hRun) runs=hRun->GetEntries();
01314   else cout<<"Ahhh, no hRun histo"<<endl;
01315 
01316   //get the spills per file
01317   Float_t spillsPerFile=-1;
01318   s="hSpillsPerFile";
01319   s+=sNamePostfix;
01320   TH1F* hSpillsPerFile=(TH1F*)gROOT->FindObjectAny(s.c_str());
01321   if (hSpillsPerFile){
01322     spillsPerFile=hSpillsPerFile->GetMean();
01323   }
01324   else {
01325     spillsPerFile=400;
01326     cout<<"Can't find hSpillsPerFile, so using spillsPerFile="
01327         <<spillsPerFile<<endl;
01328 
01329     //deal with a bug
01330     if (spillsPerFile==0 && det==Detector::kFar) {
01331       spillsPerFile=1;
01332       cout<<endl<<"*** WARNING spillsPerFile=0 so setting to "
01333           <<spillsPerFile<<" ***"<<endl;
01334     }
01335   }
01336 
01337   //apply the necessary 1e12 factor
01338   avPot*=1e12;
01339 
01340   //for the FD there is an additional factor of 1e8
01341   if (det==Detector::kFar) avPot*=1e8;
01342 
01343   Float_t totalPotFiles=0;//calc from number of files
01344   if (simFlag==SimFlag::kMC) {
01345     cout<<"For spills per files method: using avPot="<<avPot
01346         <<", avPotBad="<<avPotBad
01347         <<", runs="<<runs<<", spillsPerFile="<<spillsPerFile<<endl;
01348 
01349     //calc the totalPotFiles
01350     totalPotFiles=avPot*spillsPerFile*runs;
01351     cout<<"Total POT calc from spills per file="
01352         <<totalPotFiles<<endl;
01353   }
01354 
01355   //now calc the pots using the information from individual spills
01356   //this is for data
01357   Float_t totalPotHist=0;
01358   Float_t totalPotBadHist=0;
01359   if (det==Detector::kNear && simFlag==SimFlag::kData) {
01360     if (hPot) totalPotHist=
01361                 hPot->GetMean()*hPot->GetEntries();
01362     else cout<<"Ahhh, no hPottortgt histogram"<<endl;
01363 
01364     if (hPotBad) totalPotBadHist=
01365                    hPotBad->GetMean()*hPotBad->GetEntries();
01366     else cout<<"Ahhh, no hPotBadtortgt histogram"<<endl;
01367 
01368     //multiply by the necessary 1e12 factor
01369     totalPotHist*=1e12;
01370     totalPotBadHist*=1e12;
01371 
01372     cout<<"Total POT calc from spill beam data="<<totalPotHist
01373         <<" (potsRejected="<<totalPotBadHist<<")"<<endl;
01374   }
01375 
01376   //decide which version to use
01377   Float_t totalPot=1;
01378   if (simFlag==SimFlag::kData) {
01379     totalPot=totalPotHist;
01380   }
01381   else if (simFlag==SimFlag::kMC) {
01382     totalPot=totalPotFiles;
01383   }
01384   else cout<<"ahhh, pot error"<<endl;
01385 
01387   //section to create the hTotalPot histo (or not)
01389   if (simFlag==-1) {
01390     //do nothing
01391     cout<<"Not creating hTotalPot histo"<<endl;
01392   }
01393   else if (det==Detector::kFar && simFlag==SimFlag::kData) {
01394     //don't create the histogram here
01395     //FD data pot counting has to be done by hand
01396     cout<<"Not creating hTotalPot histo for det="<<det
01397         <<", simFlag="<<simFlag<<endl;
01398   }
01399   else {
01400     //make a histogram to store the result
01401     TH1F* hTotalPot=new TH1F("hTotalPot","hTotalPot",2,-0.5,1.5);
01402     hTotalPot->Fill(1,totalPot);
01403     cout<<"hTotalPot->GetMax()="<<hTotalPot->GetMaximum()
01404         <<", hTotalPot->Integral()="<<hTotalPot->Integral()<<endl;
01405   }
01406 
01407   return totalPot;
01408 }
01409 
01410 //......................................................................