NuUtilities Class Reference

#include <NuUtilities.h>

List of all members.

Public Member Functions

 NuUtilities ()
virtual ~NuUtilities ()
 ClassDef (NuUtilities, 0)

Static Public Member Functions

static TGraph ConvertLikelihoodSurfaceToContour (TH2D *hLikelihood, Double_t upValue)
static Double_t DistanceTargetToND ()
static void CalcFiducialMass ()
static TString PrintRelease (int releaseType)
static BeamType::BeamType_t ReturnConventionsBeamType (const NuConfig &config)
static BeamType::BeamType_t ConvertToConventionsBeamType (const Int_t beamType, const Int_t intensity)
static std::vector< Double_t > GenerateBins (int ngroups, double edges[], double steps[])
static std::vector< Double_t > RecoBins (const NuBinningScheme::NuBinningScheme_t binningScheme)
static std::vector< Double_t > TrueBins (const NuBinningScheme::NuBinningScheme_t binningScheme)
static Double_t ValueAt (TH1 *h, Double_t val)
static void WriteHistosToFile (const std::vector< TH1D > &vHistos, TFile &file, bool Quiet=false)
static void WriteHistosToFile (const std::vector< NuMatrixSpectrum > &vHistos, TFile &file, bool Quiet=false)
static void WriteHistosToFile (const std::pair< NuMatrixSpectrum, NuMatrixSpectrum > &vHistos, TFile &file)
static TH1D * RebinHistogram (const TH1D *old, int nbins, const double binedges[999])
static TH1D * RebinHistogram (const TH1D *old, std::vector< Double_t > &bins)
static TH1D * RebinHistogram (const TH1D *old, NuBinningScheme::NuBinningScheme_t scheme)
static Double_t LogLikelihood (Double_t data, Double_t MC)
static Double_t GetChisqFromCL (double cl=0.6827, int n=1)
static Double_t GetChisqFromSigma (double s=1, int n=1)
static void OscillateJesstogram (TH2D *h, double dm2, double sn2)
static void FixDogwoodQP (NuEvent &nu)
static std::vector< TString > GetListOfFilesInDir (const TString &wildcardString)
static int GetDigit (int num, int d)
static bool CheckConsistency (const TH1 *h1, const TH1 *h2)
 Checks two histograms have the same binning.
static const TH1D IntegralMarginalize (const TH2D *surf, const char axis= 'y')
static const TH1D MinimumMarginalize (const TH2D *surf, const char axis= 'y')
static void ProgressBar (Int_t entry, Float_t nEntriesF, Int_t mintime=0)

Static Public Attributes

static VldTimeStamp vtsRunIEnd
static VldTimeStamp vtsRunIIEnd
static VldTimeStamp vtsRunIIIEnd
static VldTimeStamp vtsRunIVEnd
static VldTimeStamp vtsRunVEnd
static VldTimeStamp vtsRunVIEnd
static VldTimeStamp vtsRunVIIEnd
static VldTimeStamp vtsRunVIIIEnd
static VldTimeStamp vtsRunIXEnd
static VldTimeStamp vtsRunXEnd


Detailed Description

Definition at line 46 of file NuUtilities.h.


Constructor & Destructor Documentation

NuUtilities::NuUtilities (  ) 

Definition at line 55 of file NuUtilities.cxx.

00056 {
00057 }

NuUtilities::~NuUtilities (  )  [virtual]

Definition at line 60 of file NuUtilities.cxx.

00061 {
00062 }


Member Function Documentation

void NuUtilities::CalcFiducialMass (  )  [static]

Definition at line 82 of file NuUtilities.cxx.

References Munits::cm, Munits::cm3, Munits::g, Munits::kilotonne, Munits::m, and Munits::m2.

00083 {
00084 
00085 
00086   const Float_t avPlSteelDensityNDData=7.847*(Munits::g/Munits::cm3);
00087   const Float_t avPlSteelThicknessNDData=2.563*(Munits::cm);
00088   const Float_t avPlScintDensityNDData=1.0457*(Munits::g/Munits::cm3);
00089   const Float_t avPlScintThicknessNDData=1.02*(Munits::cm);
00090   const Float_t avPlAlumDensityNDData=2.7*(Munits::g/Munits::cm3);
00091   const Float_t avPlAlumThicknessNDData=0.1*(Munits::cm);
00092 
00093   const Float_t avPlSteelDensityFDData=avPlSteelDensityNDData;
00094   const Float_t avPlSteelThicknessFDData=2.559*(Munits::cm);
00095   const Float_t avPlScintDensityFDData=avPlScintDensityNDData;
00096   const Float_t avPlScintThicknessFDData=avPlScintThicknessNDData;
00097   const Float_t avPlAlumDensityFDData=avPlAlumDensityNDData;
00098   const Float_t avPlAlumThicknessFDData=avPlAlumThicknessNDData;
00099 
00100   //MC
00101   const Float_t avPlSteelDensityNDMC=avPlSteelDensityNDData;
00102   const Float_t avPlSteelThicknessNDMC=2.54*(Munits::cm);
00103   const Float_t avPlScintDensityNDMC=avPlScintDensityNDData;
00104   const Float_t avPlScintThicknessNDMC=1.0*(Munits::cm);
00105   const Float_t avPlAlumDensityNDMC=avPlAlumDensityNDData;
00106   const Float_t avPlAlumThicknessNDMC=avPlAlumThicknessNDData;
00107 
00108   const Float_t avPlSteelDensityFDMC=avPlSteelDensityNDMC;
00109   const Float_t avPlSteelThicknessFDMC=avPlSteelThicknessNDMC;
00110   const Float_t avPlScintDensityFDMC=avPlScintDensityNDMC;
00111   const Float_t avPlScintThicknessFDMC=avPlScintThicknessNDMC;
00112   const Float_t avPlAlumDensityFDMC=avPlAlumDensityNDMC;
00113   const Float_t avPlAlumThicknessFDMC=avPlAlumThicknessNDMC;
00114 
00115   cout<<"************** Thickness and Density ******************"<<endl;
00116   cout<<"Thickness ND Data:"
00117       <<" steel="<<avPlSteelThicknessNDData
00118       <<", scint="<<avPlScintThicknessNDData
00119       <<", alum="<<avPlAlumThicknessNDData<<endl;
00120   cout<<"Density ND Data:"
00121       <<" steel="<<avPlSteelDensityNDData
00122       <<", scint="<<avPlScintDensityNDData
00123       <<", alum="<<avPlAlumDensityNDData<<endl;
00124   cout<<endl;
00125     cout<<"Thickness ND MC:"
00126       <<" steel="<<avPlSteelThicknessNDMC
00127       <<", scint="<<avPlScintThicknessNDMC
00128       <<", alum="<<avPlAlumThicknessNDMC<<endl;
00129   cout<<"Density ND MC:"
00130       <<" steel="<<avPlSteelDensityNDMC
00131       <<", scint="<<avPlScintDensityNDMC
00132       <<", alum="<<avPlAlumDensityNDMC<<endl;
00133 
00134   cout<<endl<<endl;
00135   cout<<"Thickness FD Data:"
00136       <<" steel="<<avPlSteelThicknessFDData
00137       <<", scint="<<avPlScintThicknessFDData
00138       <<", alum="<<avPlAlumThicknessFDData<<endl;
00139   cout<<"Density FD Data:"
00140       <<" steel="<<avPlSteelDensityFDData
00141   <<", scint="<<avPlScintDensityFDData
00142       <<", alum="<<avPlAlumDensityFDData<<endl;
00143   cout<<endl;
00144     cout<<"Thickness FD MC:"
00145       <<" steel="<<avPlSteelThicknessFDMC
00146       <<", scint="<<avPlScintThicknessFDMC
00147       <<", alum="<<avPlAlumThicknessFDMC<<endl;
00148   cout<<"Density FD MC:"
00149       <<" steel="<<avPlSteelDensityFDMC
00150       <<", scint="<<avPlScintDensityFDMC
00151       <<", alum="<<avPlAlumDensityFDMC<<endl;
00152 
00153   const Float_t avPlSteelStPwNDData=
00154     avPlSteelThicknessNDData*avPlSteelDensityNDData;
00155   const Float_t avPlScintStPwNDData=
00156     avPlScintThicknessNDData*avPlScintDensityNDData;
00157   const Float_t avPlAlumStPwNDData=
00158     avPlAlumThicknessNDData*avPlAlumDensityNDData;
00159 
00160   const Float_t avPlSteelStPwFDData=
00161     avPlSteelThicknessFDData*avPlSteelDensityFDData;
00162   const Float_t avPlScintStPwFDData=
00163     avPlScintThicknessFDData*avPlScintDensityFDData;
00164   const Float_t avPlAlumStPwFDData=
00165     avPlAlumThicknessFDData*avPlAlumDensityFDData;
00166 
00167   const Float_t avPlSteelStPwNDMC=
00168     avPlSteelThicknessNDMC*avPlSteelDensityNDMC;
00169   const Float_t avPlScintStPwNDMC=
00170     avPlScintThicknessNDMC*avPlScintDensityNDMC;
00171   const Float_t avPlAlumStPwNDMC=
00172     avPlAlumThicknessNDMC*avPlAlumDensityNDMC;
00173 
00174   const Float_t avPlSteelStPwFDMC=
00175     avPlSteelThicknessFDMC*avPlSteelDensityFDMC;
00176   const Float_t avPlScintStPwFDMC=
00177     avPlScintThicknessFDMC*avPlScintDensityFDMC;
00178   const Float_t avPlAlumStPwFDMC=
00179     avPlAlumThicknessFDMC*avPlAlumDensityFDMC;
00180 
00181   cout<<endl<<endl;
00182   cout<<"**************** Stopping Power ********************"<<endl;
00183   cout<<"StPw ND Data:"
00184       <<" steel="<<avPlSteelStPwNDData
00185       <<", scint="<<avPlScintStPwNDData
00186       <<", alum="<<avPlAlumStPwNDData<<endl;
00187   cout<<endl;
00188   cout<<"StPw ND MC:"
00189       <<" steel="<<avPlSteelStPwNDMC
00190       <<", scint="<<avPlScintStPwNDMC
00191       <<", alum="<<avPlAlumStPwNDMC<<endl;
00192   cout<<endl<<endl;
00193   cout<<"StPw FD Data:"
00194       <<" steel="<<avPlSteelStPwFDData
00195       <<", scint="<<avPlScintStPwFDData
00196       <<", alum="<<avPlAlumStPwFDData<<endl;
00197   cout<<endl;
00198   cout<<"StPw FD MC:"
00199       <<" steel="<<avPlSteelStPwFDMC
00200       <<", scint="<<avPlScintStPwFDMC
00201       <<", alum="<<avPlAlumStPwFDMC<<endl;
00202 
00203 
00204   //NOW CALC AREAS
00205   Float_t areaPlND=3.142*pow(0.8*Munits::m,2);
00206   Float_t areaPlFD=3.142*(14*Munits::m2-pow(0.4*Munits::m,2));
00207   cout<<endl<<endl;
00208   cout<<"Area ND = "<<areaPlND<<endl;
00209   cout<<"Area FD = "<<areaPlFD<<endl;
00210 
00211   //NOW CALC MASS
00212   const Float_t avPlSteelMassNDData=avPlSteelStPwNDData*areaPlND;
00213   const Float_t avPlScintMassNDData=avPlScintStPwNDData*areaPlND;
00214   const Float_t avPlAlumMassNDData=avPlAlumStPwNDData*areaPlND;
00215 
00216   const Float_t avPlSteelMassNDMC=avPlSteelStPwNDMC*areaPlND;
00217   const Float_t avPlScintMassNDMC=avPlScintStPwNDMC*areaPlND;
00218   const Float_t avPlAlumMassNDMC=avPlAlumStPwNDMC*areaPlND;
00219 
00220   const Float_t avPlSteelMassFDData=avPlSteelStPwFDData*areaPlFD;
00221   const Float_t avPlScintMassFDData=avPlScintStPwFDData*areaPlFD;
00222   const Float_t avPlAlumMassFDData=avPlAlumStPwFDData*areaPlFD;
00223 
00224   const Float_t avPlSteelMassFDMC=avPlSteelStPwFDMC*areaPlFD;
00225   const Float_t avPlScintMassFDMC=avPlScintStPwFDMC*areaPlFD;
00226   const Float_t avPlAlumMassFDMC=avPlAlumStPwFDMC*areaPlFD;
00227 
00228   cout<<endl<<endl;
00229   cout<<"**************** MASS ********************"<<endl;
00230   cout<<"Mass ND Data:"
00231       <<" steel="<<avPlSteelMassNDData
00232       <<", scint="<<avPlScintMassNDData
00233       <<", alum="<<avPlAlumMassNDData<<endl;
00234   cout<<endl;
00235   cout<<"Mass ND MC:"
00236       <<" steel="<<avPlSteelMassNDMC
00237       <<", scint="<<avPlScintMassNDMC
00238       <<", alum="<<avPlAlumMassNDMC<<endl;
00239   cout<<endl<<endl;
00240   cout<<"Mass FD Data:"
00241       <<" steel="<<avPlSteelMassFDData
00242       <<", scint="<<avPlScintMassFDData
00243       <<", alum="<<avPlAlumMassFDData<<endl;
00244   cout<<endl;
00245   cout<<"Mass FD MC:"
00246       <<" steel="<<avPlSteelMassFDMC
00247       <<", scint="<<avPlScintMassFDMC
00248       <<", alum="<<avPlAlumMassFDMC<<endl;
00249 
00250 
00251 
00252   const Float_t avPlTotalMassNDData=
00253     avPlSteelMassNDData+avPlScintMassNDData+avPlAlumMassNDData;
00254   const Float_t avPlTotalMassNDMC=
00255     avPlSteelMassNDMC+avPlScintMassNDMC+avPlAlumMassNDMC;
00256   const Float_t avPlTotalMassFDData=
00257     avPlSteelMassFDData+avPlScintMassFDData+avPlAlumMassFDData;
00258   const Float_t avPlTotalMassFDMC=
00259     avPlSteelMassFDMC+avPlScintMassFDMC+avPlAlumMassFDMC;
00260 
00261   const Float_t massSumND=avPlTotalMassNDData+avPlTotalMassNDMC;
00262   const Float_t massSumFD=avPlTotalMassFDData+avPlTotalMassFDMC;
00263   const Float_t massDiffND=avPlTotalMassNDData-avPlTotalMassNDMC;
00264   const Float_t massDiffFD=avPlTotalMassFDData-avPlTotalMassFDMC;
00265 
00266   const Float_t massRelDiffND=massDiffND/(0.5*massSumND);
00267   const Float_t massRelDiffFD=massDiffFD/(0.5*massSumFD);
00268 
00269   cout<<endl<<endl;
00270   cout<<"Mass ND Data Total: "<<avPlTotalMassNDData<<endl;
00271   cout<<"Mass ND MC Total: "<<avPlTotalMassNDMC<<endl;
00272   cout<<"Mass ND Data-MC="<<massDiffND<<endl;
00273   cout<<"Mass ND Data-MC/(0.5*Data+MC)="<<massRelDiffND<<endl;
00274   cout<<endl;
00275   cout<<"Mass FD Data Total: "<<avPlTotalMassFDData<<endl;
00276   cout<<"Mass FD MC Total: "<<avPlTotalMassFDMC<<endl;
00277   cout<<"Mass FD Data-MC="<<massDiffFD<<endl;
00278   cout<<"Mass FD Data-MC/(0.5*Data+MC)="<<massRelDiffFD<<endl;
00279   cout<<endl;
00280 
00281   //Considering events from STEEL planes
00282   //ND: 14->68 inclusive = 55 planes of steel
00283   cout<<endl<<endl;
00284   cout<<"55 planes in ND Data="
00285       <<55*avPlTotalMassNDData/(Munits::kilotonne)<<" kilotonnes"<<endl;
00286   cout<<"55 planes in ND MC  ="
00287       <<55*avPlTotalMassNDMC/(Munits::kilotonne)<<" kilotonnes"<<endl;
00288 
00289   //FD: 4->239 inclusive + 253->464 inclusive = 236 + 212 = 448 planes
00290   cout<<endl<<endl;
00291   cout<<"448 planes in FD Data="
00292       <<448*avPlTotalMassFDData/(Munits::kilotonne)<<" kilotonnes"
00293       <<endl;
00294   cout<<"448 planes in FD MC  ="
00295       <<448*avPlTotalMassFDMC/(Munits::kilotonne)<<" kilotonnes"
00296       <<endl;
00297 
00298 
00299 //Table 5 of minos-doc-3134v2 has updated numbers
00300 //Both detectors have the same densities in MC
00301 //steel = 7.847
00302 //scint = 1.0457, thickness=1.0cm
00303 //The pitches of ND and FD are slightly different
00304 
00305 //Peter Litchfield's email of 20/06/07 22:27 GMT
00306 //Real detectors
00307 //Plane pitch  FD (average) 5.946cm,  ND 5.94cm (3134)
00308 //Steel density  Both Detectors 7.847gm/cc (1431)
00309 //Steel thickness FD 2.559cm, ND 2.563cm (1529)
00310 //Scintillator 1.02cm density 1.0457gm/cc (scint+coating, 3134)
00311 //Aluminium    0.1cm  density 2.7gm/cc (PDG book)
00312 //FD average density = 3.602gm/cc
00313 //ND average density = 3.611gm/cc
00314 
00315 //In the Birch MC (3134)
00316 //Plane pitch  FD (average) 5.946cm,  ND 5.94cm
00317 //Steel density 7.87gm/cc
00318 //Steel thickness 2.54cm
00319 //Scintillator 1.02cm density 1.032gm/cc
00320 //Aluminium 0.1cm air, negligible density
00321 //FD avearge density = 3.539gm/cc
00322 //ND average density = 3.542gm/cc
00323 }

bool NuUtilities::CheckConsistency ( const TH1 *  h1,
const TH1 *  h2 
) [static]

Checks two histograms have the same binning.

This is a static function on TH1, but for some reason it's private Copied the implemention here (and made it saner) as the lesser of two evils.

(The greater of two evils, preserved for posterity, is as follows)

define CheckConsistency fgHackVar; public: static bool CheckConsistency include "TH1D.h" undef CheckConsistency

Definition at line 1406 of file NuUtilities.cxx.

References CheckBinLimits().

Referenced by NuMatrixSpectrum::FillOscCache(), and NuMatrix1D::FillOscCache().

01407 {
01408   if(h1->GetNbinsX() != h2->GetNbinsX() ||
01409      h1->GetNbinsY() != h2->GetNbinsY() ||
01410      h1->GetNbinsZ() != h2->GetNbinsZ())
01411     return false;
01412 
01413   if(h1->GetXaxis()->GetXmin() != h2->GetXaxis()->GetXmin() ||
01414      h1->GetXaxis()->GetXmax() != h2->GetXaxis()->GetXmax() ||
01415      h1->GetYaxis()->GetXmin() != h2->GetYaxis()->GetXmin() ||
01416      h1->GetYaxis()->GetXmax() != h2->GetYaxis()->GetXmax() ||
01417      h1->GetZaxis()->GetXmin() != h2->GetZaxis()->GetXmin() ||
01418      h1->GetZaxis()->GetXmax() != h2->GetZaxis()->GetXmax())
01419     return false;
01420 
01421   if(!CheckBinLimits(h1->GetXaxis()->GetXbins(), h2->GetXaxis()->GetXbins()) ||
01422      !CheckBinLimits(h1->GetYaxis()->GetXbins(), h2->GetYaxis()->GetXbins()) ||
01423      !CheckBinLimits(h1->GetZaxis()->GetXbins(), h2->GetZaxis()->GetXbins()))
01424     return false;
01425 
01426   return true;
01427 }

NuUtilities::ClassDef ( NuUtilities  ,
 
)

TGraph NuUtilities::ConvertLikelihoodSurfaceToContour ( TH2D *  hLikelihood,
Double_t  upValue 
) [static]

Definition at line 65 of file NuUtilities.cxx.

00067 {
00068   Double_t minimumLike = hLikelihood->GetMinimum();
00069   hLikelihood->SetContour(1);
00070   hLikelihood->SetContourLevel(0, minimumLike + upValue);
00071   TCanvas canvas;
00072   hLikelihood->Draw("contlist same");
00073   canvas.Update(); //This line of code embodies everything I hate about Root
00074 
00075   TObjArray* contArr = (TObjArray*) gROOT->GetListOfSpecials()->FindObject("contours"); //What was Rene Brun thinking?
00076   TList* contList = (TList*) contArr->At(0);
00077   //  TGraph* contGraph68 = new TGraph(*((TGraph*)contList->First()));
00078   return *((TGraph*)contList->First());
00079 }

BeamType::BeamType_t NuUtilities::ConvertToConventionsBeamType ( const Int_t  beamType,
const Int_t  intensity 
) [static]

Definition at line 360 of file NuUtilities.cxx.

References BeamType::kL010z000i, BeamType::kL010z000i_i209, BeamType::kL010z000i_i225, BeamType::kL010z000i_i232, BeamType::kL010z000i_i259, BeamType::kL010z000i_i300, BeamType::kL010z000i_i317, BeamType::kL010z000i_i326, BeamType::kL010z000i_i380, BeamType::kL010z185i, BeamType::kL010z185i_i124, BeamType::kL010z185i_i191, BeamType::kL010z185i_i213, BeamType::kL010z185i_i224, BeamType::kL010z185i_i232, BeamType::kL010z185i_i243, BeamType::kL010z185i_i257, BeamType::kL010z185i_i282, BeamType::kL010z185i_i303, BeamType::kL010z185i_i324, BeamType::kL010z185i_i361, BeamType::kL010z185i_rev, BeamType::kL010z185i_rev_i277, BeamType::kL010z185i_rev_i306, BeamType::kL010z185i_rev_i312, BeamType::kL010z185i_rev_i314, BeamType::kL010z185i_rev_i348, BeamType::kL010z185i_rev_i349, BeamType::kL010z185i_rev_i353, BeamType::kL250z200i, BeamType::kL250z200i_i100, BeamType::kL250z200i_i114, BeamType::kL250z200i_i130, BeamType::kL250z200i_i152, BeamType::kL250z200i_i165, BeamType::kL250z200i_i194, BeamType::kL250z200i_i232, BeamType::kUnknown, Msg::kWarning, and MSG.

Referenced by ReturnConventionsBeamType().

00361 {
00362   BeamType::BeamType_t originalBeamType =
00363     static_cast<BeamType::BeamType_t>(beamType);
00364   if (0==intensity){
00365     return originalBeamType;
00366   }
00367   if (BeamType::kL010z185i == originalBeamType){
00368     if (124==intensity){
00369       return BeamType::kL010z185i_i124;
00370     }
00371     else if (191==intensity){
00372       return BeamType::kL010z185i_i191;
00373     }
00374     else if (213==intensity){
00375       return BeamType::kL010z185i_i213;
00376     }
00377     else if (224==intensity){
00378       return BeamType::kL010z185i_i224;
00379     }
00380     else if (232==intensity){
00381       return BeamType::kL010z185i_i232;
00382     }
00383     else if (243==intensity){
00384       return BeamType::kL010z185i_i243;
00385     }
00386     else if (257==intensity){
00387       return BeamType::kL010z185i_i257;
00388     }
00389     else if (282==intensity){
00390       return BeamType::kL010z185i_i282;
00391     }
00392     else if (303==intensity){
00393       return BeamType::kL010z185i_i303;
00394     }
00395     else if (324==intensity){
00396       return BeamType::kL010z185i_i324;
00397     }
00398     else if (361==intensity){
00399       return BeamType::kL010z185i_i361;
00400     }
00401     else {
00402       MSG("NuUtilities",Msg::kWarning)
00403         << "Received an L010z185i beam type with unrecognized intensity. Intensity = "
00404         << intensity
00405         << endl;
00406       return BeamType::kUnknown;
00407     }
00408   }
00409 
00410   if (BeamType::kL010z000i == originalBeamType){
00411     if (209 == intensity){
00412       return BeamType::kL010z000i_i209;
00413     }
00414     else if (225 == intensity){
00415       return BeamType::kL010z000i_i225;
00416     }
00417     else if (232 == intensity){
00418       return BeamType::kL010z000i_i232;
00419     }
00420     else if (259 == intensity){
00421       return BeamType::kL010z000i_i259;
00422     }
00423     else if (300 == intensity){
00424       return BeamType::kL010z000i_i300;
00425     }
00426     else if (317 == intensity){
00427       return BeamType::kL010z000i_i317;
00428     }
00429     else if (326 == intensity){
00430       return BeamType::kL010z000i_i326;
00431     }
00432     else if (380 == intensity){
00433       return BeamType::kL010z000i_i380;
00434     }
00435     else {
00436       MSG("NuUtilities",Msg::kWarning)
00437         << "Received an L010z000i beam type with unrecognized intensity."
00438         << endl;
00439       return BeamType::kUnknown;
00440     }
00441   }
00442 
00443   if (BeamType::kL250z200i == originalBeamType){
00444     if (100 == intensity){
00445       return BeamType::kL250z200i_i100;
00446     }
00447     else if (114 == intensity){
00448       return BeamType::kL250z200i_i114;
00449     }
00450     else if (130 == intensity){
00451       return BeamType::kL250z200i_i130;
00452     }
00453     else if (152 == intensity){
00454       return BeamType::kL250z200i_i152;
00455     }
00456     else if (165 == intensity){
00457       return BeamType::kL250z200i_i165;
00458     }
00459     else if (194 == intensity){
00460       return BeamType::kL250z200i_i194;
00461     }
00462     else if (232 == intensity){
00463       return BeamType::kL250z200i_i232;
00464     }
00465     else {
00466       MSG("NuUtilities",Msg::kWarning)
00467         << "Received an L250z200i beam type with unrecognized intensity."
00468         << endl;
00469       return BeamType::kUnknown;
00470     }
00471   }
00472 
00473   if(BeamType::kL010z185i_rev == originalBeamType){
00474     if(277 == intensity){
00475       return BeamType::kL010z185i_rev_i277;
00476     }
00477     else if(314 == intensity){
00478       return BeamType::kL010z185i_rev_i314;
00479     }
00480     else if(349 == intensity){
00481       return BeamType::kL010z185i_rev_i349;
00482     }
00483     else if(312 == intensity){
00484       return BeamType::kL010z185i_rev_i312;
00485     }
00486     else if(348 == intensity){
00487       return BeamType::kL010z185i_rev_i348;
00488     }
00489     else if(306 == intensity){
00490       return BeamType::kL010z185i_rev_i306;
00491     }
00492     else if(353 == intensity){
00493       return BeamType::kL010z185i_rev_i353;
00494     }
00495     else {
00496       MSG("NuUtilities",Msg::kWarning)
00497         << "Received an L010z185i_rev beam type with unrecognized intensity."
00498         << endl;
00499       return BeamType::kUnknown;
00500     }
00501   }
00502 
00503   MSG("NuUtilities",Msg::kWarning)
00504     << "Received a non-zero intensity flag, along with a beam "
00505     << "type not configured for intensity weighting."
00506     << endl;
00507   return BeamType::kUnknown;
00508 
00509 }

static Double_t NuUtilities::DistanceTargetToND (  )  [inline, static]

Definition at line 54 of file NuUtilities.h.

References Munits::km.

Referenced by NuOscProb::FourFlavourDisappearanceWeight(), NuOscProb::FourFlavourNuESurvivalProbability(), NuOscProb::FourFlavourNuMuToNuEProbability(), NuOscProb::FourFlavourNuMuToNuSProbability(), NuOscProb::FourFlavourNuMuToNuTauProbability(), and NuMatrix1D::OscillatedNDRecoToTrue().

00054 {return 1.04*Munits::km;};

void NuUtilities::FixDogwoodQP ( NuEvent nu  )  [static]

Definition at line 1239 of file NuUtilities.cxx.

References NuReco::GetEvtEnergy(), NuLibrary::Instance(), ReleaseType::IsDogwood(), Msg::kInfo, MAXMSG, PrintRelease(), NuEvent::qp, NuEvent::qp_rangebiased, NuEvent::qp_sigqp, NuLibrary::reco, NuEvent::recoVersion, NuEvent::sigqp, and NuEvent::sigqp_qp.

Referenced by NuDSTAna::QPStudy(), and NuDSTAna::SelectorTable().

01240 {
01241   // Only need to fix dogwood NuEvent's
01242   if (!ReleaseType::IsDogwood(nu.recoVersion)) {
01243     MAXMSG("NuUtilities",Msg::kInfo,5) << "Not correcting q/p for "
01244     << PrintRelease(nu.recoVersion) << endl;
01245     return;
01246   }
01247   MAXMSG("NuUtilities",Msg::kInfo,5) << "Correcting "
01248   << PrintRelease(nu.recoVersion) << " q/p from "
01249   << nu.qp << " to " << nu.qp_rangebiased << endl;
01250 
01251   nu.qp = nu.qp_rangebiased;
01252   if (nu.sigqp == 0) nu.qp_sigqp = 0;
01253   else               nu.qp_sigqp = nu.qp_rangebiased / nu.sigqp;
01254 
01255   for(int trkIdx = 1; trkIdx <= 3; ++trkIdx){
01256     get_qp(nu, trkIdx) = get_qp_rangebiased(nu, trkIdx);
01257     if(get_sigqp(nu, trkIdx) != 0){
01258       get_qp_sigqp(nu, trkIdx) =
01259         get_qp_rangebiased(nu, trkIdx) / get_sigqp(nu, trkIdx);
01260     }
01261   }
01262 
01263   if (nu.qp_sigqp == 0) nu.sigqp_qp = 0;
01264   else                  nu.sigqp_qp = 1./nu.qp_sigqp;
01265 
01266   NuLibrary& lib=NuLibrary::Instance();
01267   lib.reco.GetEvtEnergy(nu, false);
01268 
01269   if (nu.qp_rangebiased > 0) nu.charge = 1;
01270   else                       nu.charge = -1;
01271 }

vector< Double_t > NuUtilities::GenerateBins ( int  ngroups,
double  edges[],
double  steps[] 
) [static]

Definition at line 512 of file NuUtilities.cxx.

References Munits::g.

Referenced by NuMatrixSpectrum::RebinForFit(), NuMatrixSpectrum::RebinForPlotAsFit(), NuMatrixSpectrum::RebinForPlots(), NuMatrixSpectrum::RebinToGeV(), and RecoBins().

00513 {
00514     vector<Double_t> bins;
00515     for (int g = 0; g < ngroups; g++) {
00516         Double_t lowEdge = edges[g];
00517         Double_t highEdge = edges[g+1];
00518         Double_t binWidth = steps[g];
00519         for (Double_t binEdge = lowEdge;
00520              binEdge < highEdge - (binWidth/2.0);
00521              binEdge += binWidth){
00522             bins.push_back(binEdge);
00523         }
00524     }
00525     bins.push_back(edges[ngroups]);
00526     return bins;
00527 }

Double_t NuUtilities::GetChisqFromCL ( double  cl = 0.6827,
int  n = 1 
) [static]

Definition at line 1337 of file NuUtilities.cxx.

Referenced by GetChisqFromSigma().

01338 {
01339 
01340   if(n<2) return 2*pow(TMath::ErfInverse(cl),2);
01341 
01342    // Create the function and wrap it
01343    TF1 f("C.L. from Chisq", "TMath::Gamma(0.5*[0],0.5*x)", 0, 2*n+50);
01344    f.SetParameter(0,n);
01345    f.SetNpx(4);
01346  
01347    return f.GetX(cl);
01348 
01349 }

Double_t NuUtilities::GetChisqFromSigma ( double  s = 1,
int  n = 1 
) [static]

Definition at line 1353 of file NuUtilities.cxx.

References GetChisqFromCL(), and Munits::s.

01354 {
01355 
01356   if(n<2) return s*s;
01357 
01358   double cl = TMath::Erf(s/sqrt(2));
01359   
01360   return GetChisqFromCL(cl,n);
01361 
01362 }

int NuUtilities::GetDigit ( int  num,
int  d 
) [static]

Definition at line 1384 of file NuUtilities.cxx.

Referenced by NuExtraction::ExtractCoilInfo().

01385 {
01386   int pow1 = static_cast<int>(TMath::Power(10, d));
01387   int pow2 = pow1 / 10;
01388   return (num % pow1 - num % pow2) / pow2;
01389 }

std::vector< TString > NuUtilities::GetListOfFilesInDir ( const TString &  wildcardString  )  [static]

Definition at line 1275 of file NuUtilities.cxx.

References Munits::g.

Referenced by NuFluxHelper::ConcatenateHelpers(), and NuFluxHelper::MakeHelperHistos().

01276 {
01277   std::vector<TString> fileList;
01278 
01279   glob_t g;
01280   glob(wildcardString.Data(),
01281        // GLOB_NOCHECK | // If no match return pattern
01282        GLOB_TILDE, // Expand ~'s
01283        0, &g);
01284   for(unsigned int i = 0; i < g.gl_pathc; ++i)
01285     fileList.push_back(g.gl_pathv[i]);
01286   globfree(&g);
01287 
01288   // glob output is sorted by default
01289 
01290   return fileList;
01291 }

const TH1D NuUtilities::IntegralMarginalize ( const TH2D *  surf,
const char  axis = 'y' 
) [static]

Definition at line 1430 of file NuUtilities.cxx.

References Msg::kError, and MSG.

01431 {
01432   bool doY;
01433   if (axis == 'y' || axis == 'Y')  doY = true;
01434   else if (axis == 'x' || axis == 'X') doY = false;
01435   else {
01436     MSG("NuUtilities",Msg::kError) << "Axis " << axis << "not recognized. Only x, X, y, Y are allowed axes." << endl;
01437     TH1D blank;    return blank;
01438   }
01439 
01440   TH1D * marginSurf;
01441   TString name = surf->GetName();
01442 
01443   // Use root's built in projections to do the integrals, excluding Under/overflow
01444   if (doY) {
01445     marginSurf = surf->ProjectionY(name+"_intMarginY", 1, surf->GetNbinsX());
01446     marginSurf->Scale(1./surf->GetNbinsX());
01447   }
01448   else {
01449     marginSurf = surf->ProjectionX(name+"_intMarginX", 1, surf->GetNbinsY());
01450     marginSurf->Scale(1./surf->GetNbinsY());
01451   }
01452 
01453   // Subtract off the minimum bin to get a deltaChisquared
01454   // Not *exactly* correct but should be pretty close in bins are
01455   // relatively small
01456   double allMin = marginSurf->GetMinimum();
01457   for (int b = 0; b <= marginSurf->GetNbinsX()+1; b++) {
01458     double val = marginSurf->GetBinContent(b) - allMin;
01459     marginSurf->SetBinContent(b, val);
01460   }
01461 
01462   return *marginSurf;
01463 }

Double_t NuUtilities::LogLikelihood ( Double_t  data,
Double_t  MC 
) [static]

Definition at line 1294 of file NuUtilities.cxx.

References MuELoss::e.

Referenced by NuMMRun::StatsLikelihood(), NuMatrix2D::StatsLikelihood(), and NuMatrix1D::StatsLikelihood().

01295 {
01296 
01297   if( isnan(mnu) || isnan(dnu) || dnu < 0 ){
01298     cout << "Error: MC or Data are not valid numbers!!!" << endl;
01299     return 1e12;
01300   }
01301 
01302   const Double_t MCError = 1e-12;
01303   if (dnu) {
01304     // Protection against cases where there is data but no MC?
01305     // This should not normally be possible, other than highly contrived
01306     // unphysical oscillation values.
01307     // Set the value to something lower than POT_Data / POT_MC, which at the
01308     // moment is usually 1e-3. We can justify this as an error on the MC
01309     // Set to a 2nd order expansion around MCError.
01310     // Also protects against negative MC events
01311 
01312     if (mnu < MCError){
01313       //cout << "MC Error: MC = " << mnu << endl;
01314       return 2*(mnu-dnu+dnu*log(dnu/MCError)) + dnu*(3-4*mnu/MCError+pow(mnu/MCError,2));
01315     }
01316 
01317     // We have both values. Do the proper lnL!
01318     return 2*(mnu - dnu + dnu*log(dnu/mnu));
01319   }
01320   else {
01321     // dnu is zero, mnu may or may not be zero
01322     // Protect against negative numbers by setting to a
01323     // quadratic function that matches Loglikelihood up
01324     // to 1st order at MCError and has minimum at mnu=0
01325 
01326     if (mnu < MCError){
01327       //cout << "MC Error: MC = " << mnu << endl;
01328       return mnu*mnu/MCError + MCError;
01329     }
01330 
01331     return 2*mnu;
01332   }
01333 }

const TH1D NuUtilities::MinimumMarginalize ( const TH2D *  surf,
const char  axis = 'y' 
) [static]

Definition at line 1467 of file NuUtilities.cxx.

References Msg::kError, min, and MSG.

01468 {
01469   bool doY;
01470   if (axis == 'y' || axis == 'Y')  doY = true;
01471   else if (axis == 'x' || axis == 'X') doY = false;
01472   else {
01473     MSG("NuUtilities",Msg::kError) << "Axis " << axis << "not recognized. Only x, X, y, Y are allowed axes." << endl;
01474     TH1D blank;    return blank;
01475   }
01476 
01477   TH1D * marginSurf;
01478   TString name = surf->GetName();
01479 
01480   int bKeepMax;
01481   int bMargMax;
01482   // Use root's built in projections to define the axis ranges, reset contents
01483   if (doY) {
01484     marginSurf = surf->ProjectionY(name+"_minMarginY", 1, surf->GetNbinsX());
01485     bKeepMax = surf->GetNbinsY(); // Keep Y
01486     bMargMax = surf->GetNbinsX(); // Marginalize X
01487   }
01488   else {
01489     marginSurf = surf->ProjectionX(name+"_minMarginX", 1, surf->GetNbinsY());
01490     bKeepMax = surf->GetNbinsX(); // Keep X
01491     bMargMax = surf->GetNbinsY(); // Marginalize Y
01492   }
01493   marginSurf->Reset();
01494 
01495   double allMin = 1.e10;
01496   double rowMin;
01497   int bx, by;
01498   for (int bk = 1; bk < bKeepMax; bk++) {
01499     rowMin = 1.e10;
01500     for (int bm = 1; bm < bMargMax; bm++) {
01501       if (doY) {
01502         bx = bm; // Marginalize X
01503         by = bk; // Keep Y
01504       }
01505       else {
01506         bx = bk; // Keep X
01507         by = bm; // Marginalize Y
01508       }
01509       rowMin = min(surf->GetBinContent(bx, by), rowMin);
01510       allMin = min(rowMin, allMin);
01511     }
01512     marginSurf->SetBinContent(bk, rowMin);
01513   }
01514 
01515   // Subtract off the minimum bin to get a deltaChisquared
01516   // Not *exactly* correct but should be pretty close in bins are
01517   // relatively small
01518   for (int bk = 1; bk < bKeepMax; bk++) {
01519     double val = marginSurf->GetBinContent(bk) - allMin;
01520     marginSurf->SetBinContent(bk, val);
01521   }
01522 
01523   return *marginSurf;
01524 }

void NuUtilities::OscillateJesstogram ( TH2D *  h,
double  dm2,
double  sn2 
) [static]

Definition at line 1365 of file NuUtilities.cxx.

References one(), NuMatrixSpectrum::Oscillate(), and NuMatrixSpectrum::Spectrum().

01366 {
01367   const int Nx = h->GetNbinsX();
01368   const int Ny = h->GetNbinsY();
01369 
01370   TH1D one("", "", Nx, h->GetXaxis()->GetXbins()->GetArray());
01371   for(int n = 0; n < Nx+2; ++n) one.SetBinContent(n, 1);
01372 
01373   NuMatrixSpectrum osc(one);
01374   osc.Oscillate(dm2, sn2);
01375 
01376   for(int x = 0; x < Nx+2; ++x){
01377     for(int y = 0; y < Ny+2; ++y){
01378       h->SetBinContent(x, y, h->GetBinContent(x, y)*osc.Spectrum()->GetBinContent(x));
01379     }
01380   }
01381 }

TString NuUtilities::PrintRelease ( int  releaseType  )  [static]

Definition at line 326 of file NuUtilities.cxx.

References ReleaseType::AsString().

Referenced by NuAnalysis::ExtractConfig(), FixDogwoodQP(), NuPIDInterface::GetFileNameAbID(), NuPIDInterface::GetFileNamekNNID(), NuZBeamReweight::GetWeightHelium(), NuAnalysis::LIRejectionTest(), NuDSTAna::MakeMicroDstHe(), NuPIDInterface::PrintReleaseTypeEtc(), and NuAnalysis::SetAnaFlags().

00327 {
00328   char hexnum[20];
00329   sprintf(hexnum, "%x", releaseType);
00330   TString ret = ReleaseType::AsString(releaseType);
00331   ret += " (0x";
00332   ret += hexnum;
00333   ret += ", ";
00334   ret += releaseType;
00335   ret += ")";
00336   return ret;
00337 }

void NuUtilities::ProgressBar ( Int_t  entry,
Float_t  nEntriesF,
Int_t  mintime = 0 
) [static]

Prints a progress bar with ETA. Not fully finished yet, still has issues with printing when the output is directed somewhere other than the console

Definition at line 1660 of file NuUtilities.cxx.

References Munits::bar, MuELoss::e, Msg::kInfo, and MSG.

Referenced by NuFCEventManager::AddFile(), NuSystFitter::Chi2ValleyDm2(), NuSystFitter::Chi2ValleySn2(), NuSystFitter::FillAlphaBeta(), NuSystFitter::FillCPT(), NuSystFitter::FillLikelihoodSurfaceCPT(), NuSystFitter::FillLikelihoodSurfaceGeneric(), NuSystFitter::FillLikelihoodSurfaceNSIDm2Eps(), NuSystFitter::FillLikelihoodSurfaceNSIDm2Sn2(), NuSystFitter::FillLikelihoodSurfaceNSISn2Eps(), NuSystFitter::FillLikelihoodSurfaceTransDm2BarGMu(), NuSystFitter::FillLikelihoodSurfaceTransDm2BarSn2Bar(), NuSystFitter::FillLikelihoodSurfaceTransDm2GMu(), NuSystFitter::FillLikelihoodSurfaceTransDm2GMuBoth(), NuSystFitter::FillLikelihoodSurfaceTransDm2Sn2(), NuSystFitter::FillLikelihoodSurfaceTransDm2Sn2Both(), NuSystFitter::FillLikelihoodSurfaceTransSn2BarGMu(), NuSystFitter::FillLikelihoodSurfaceTransSn2GMu(), NuSystFitter::FillLikelihoodSurfaceTransSn2GMuBoth(), ToFPlotter::MakeDataPlots(), FoverNHistos::MakeDataPlots(), DataMCPlots::MakeDataPlots(), TemplateAnalysisClass::MakePlots(), FitTree::MakeVectors(), NuDSTAna::MMRereco(), and NuDSTAna::MMTransSME().

01661 {
01662   Int_t nEntries = (Int_t)nEntriesF;
01663 
01664   // Are we streaming to a file? If so, show the old style
01665   struct stat buf;
01666   fstat(fileno(stdout), &buf);
01667   // Check if we are a file, or a pipe (i.e. in case the user tees)
01668   const bool isFile = buf.st_mode & (S_IFREG | S_IFIFO) ;
01669   if (isFile) {
01670     Float_t fract = ceil(nEntries/20.);
01671     if (ceil(((Float_t)entry)/fract)==((Float_t)entry)/fract)
01672     {
01673         MSG("NuUtilities",Msg::kInfo)
01674           <<"Fraction of loop complete: "<<entry
01675           <<"/"<<nEntries<<"  ("
01676           <<(Int_t)(100.*entry/nEntries)<<"%)"<<endl;
01677     }
01678     return;
01679   }
01680 
01681 
01682   // How wide should we make the progress bar?
01683   const Int_t width = 70;
01684   // How long is the string for entries?
01685   static Int_t countlen = -1;
01686   // How long is our progress bar?
01687   static Int_t barlen = -1;
01688   // The entry number of the next bar entry
01689   static Int_t nextbar = -1;
01690   // When did we start?
01691   static time_t starttime = 0;
01692   // when do we next update?
01693   static time_t nextupdate = 0;
01694   // If we are processing the first entry, reset everything
01695   if (entry <= 1)
01696   {
01697     // Get the new length of the entry string
01698     countlen = (Int_t)TMath::Ceil(TMath::Log10(nEntries)) + 1;
01699     nextbar = -1;
01700     starttime = time(NULL);
01701 
01702     barlen = width - 14 - countlen*2 - 1;
01703 
01704     // Don't update until we get to the minimum time
01705     nextupdate = starttime + mintime;
01706   }
01707 
01708   // Check here to see if we should update; otherwise, return
01709   // Check to see if the bar would update
01710   // or, alternatively, it is time to refresh.
01711   if ((time(NULL) < nextupdate) && (entry < nextbar || (nextbar == -1)))return;
01712   nextupdate = time(NULL) + 10;
01713 
01714   // Because this is used in several places, make it here
01715   Float_t frac = (Float_t)entry / (Float_t)nEntries;
01716 
01717   // Prepare the progress bar string
01718   TString bar;
01719   if (entry <= nEntries)
01720   {
01721     // Work out how many characters we are in
01722     Int_t numeq = TMath::FloorNint(frac*barlen);
01723 
01724     // Work out when the next bar will occur
01725     nextbar = (Int_t)((Float_t)(numeq+1)/(Float_t)barlen*nEntries);
01726     //cout << "Next bar at: " << nextbar << "        " << endl;
01727     bar = TString('=', numeq);
01728     bar += TString(' ', barlen - numeq);
01729   } else if (entry > nEntries) {
01730     // We have gone over. Oh no!
01731     bar = TString('+', barlen);
01732   } else if (entry < 0) {
01733     // Somehow, we are below zero. Handle it nonetheless
01734     bar = TString('-', barlen);
01735   }
01736 
01737 
01738   // Prepare the ETA
01739   Float_t elapsed_time = (Float_t)(time(NULL) - starttime);
01740   Float_t time_left = -60;
01741   if (frac > 1e-6) {
01742     time_left = (elapsed_time / frac) - elapsed_time;
01743   }
01744   Int_t mins, seconds;
01745   mins    = (Int_t)TMath::Floor(time_left / 60.0f);
01746   seconds = (Int_t)TMath::Floor(time_left - (Float_t)(mins*60.0f));
01747   // TString ETA;
01748   ostringstream ETA;
01749 
01750   ETA << "ETA ";
01751   if ((mins < 0 || seconds < 0) || (mins == 0 && seconds == 0)) {
01752     ETA << "--:--";
01753   } else {
01754     ETA << setfill('0') << setw(2) << mins << ":" << setw(2) << seconds;
01755   }
01756 
01757   cout << " Progress: [" << bar << "] "
01758        << setw(countlen) << entry << "/" << nEntries
01759        << " " << ETA.str()
01760        << '\r'
01761        << flush;
01762   // Move to the next line, if this is the final entry!
01763   if (entry == nEntries) {
01764     cout << endl;
01765   }
01766 }

TH1D * NuUtilities::RebinHistogram ( const TH1D *  old,
NuBinningScheme::NuBinningScheme_t  scheme 
) [static]

Definition at line 1232 of file NuUtilities.cxx.

References RebinHistogram(), and RecoBins().

01233 {
01234   vector<Double_t> bins = RecoBins(scheme);
01235   return RebinHistogram(old, bins);
01236 }

TH1D * NuUtilities::RebinHistogram ( const TH1D *  old,
std::vector< Double_t > &  bins 
) [static]

Definition at line 1226 of file NuUtilities.cxx.

References RebinHistogram().

01227 {
01228   return RebinHistogram(old, bins.size()-1, &(bins[0]));
01229 }

TH1D * NuUtilities::RebinHistogram ( const TH1D *  old,
int  nbins,
const double  binedges[999] 
) [static]

Definition at line 1164 of file NuUtilities.cxx.

References count, and MuELoss::e.

Referenced by RebinHistogram(), NuMatrixSpectrum::RebinToArray(), and NuMatrixSpectrum::RebinToScheme().

01165 {
01166   static int count = 0;
01167   TString scount = "";
01168   scount += count++;
01169 
01170   TH1D* newh = new TH1D(old->GetName() + TString("_rebin") + scount, 
01171                         old->GetTitle(), nbins, binedges);
01172   bool ouwarn = true;
01173 
01174   // Floating point comparisons of bin edges don't go well otherwise
01175   const double epsilon = 1e-10;
01176 
01177   // Now, loop over the old histogram and fill the new one
01178   for(Int_t bin = 1; bin <= old->GetNbinsX(); bin++) {
01179 
01180     double oldCent = old->GetBinCenter(bin);
01181     int newbin = newh->GetXaxis()->FindFixBin(oldCent);
01182 
01183     // Make sure that bins do not cross boundaries.
01184     if (newbin > newh->GetNbinsX() || newbin == 0) {
01185       if (ouwarn) {
01186         cout << "Old bin: " << old->GetBinLowEdge(bin) << " - " << old->GetBinLowEdge(bin+1)
01187         << " and maybe others have become under/overflow" << endl;
01188         ouwarn = false;
01189       }
01190     } else if ( old->GetBinLowEdge(bin) + epsilon < newh->GetBinLowEdge(newbin) ||
01191                 old->GetBinLowEdge(bin + 1) > newh->GetBinLowEdge(newbin + 1) + epsilon) {
01192       cout << "Bin edge conflict in bin no's " << bin << ", " << newbin 
01193            << ".  No rebinning performed."
01194            << "  old bin: " << old->GetBinLowEdge(bin) << " - " << old->GetBinLowEdge(bin + 1)
01195            << ", new bin: " << newh->GetBinLowEdge(newbin) << " - " 
01196            << newh->GetBinLowEdge(newbin + 1)
01197            << endl;
01198       TH1D* old_ret = new TH1D(*old);
01199       return old_ret;
01200     }
01201 
01202     double oldCon = old->GetBinContent(bin);
01203     double oldErr = old->GetBinError(bin);
01204     double newCon = newh->GetBinContent(newbin);
01205     double newErr = newh->GetBinError(newbin);
01206 
01207     newCon += oldCon;
01208     newErr = sqrt(newErr*newErr + oldErr*oldErr);
01209     newh->SetBinContent(newbin, newCon);
01210     newh->SetBinError(newbin, newErr);
01211   }
01212 
01213   // Validate
01214   double oldint = old->Integral();
01215   double newint = newh->Integral() + newh->GetBinContent(0)
01216                  + newh->GetBinContent(newh->GetNbinsX()+1);
01217   if (abs(oldint - newint) > 1e-5) {
01218     cout << "Rebinning report: " << oldint << " entries have become " << newint << " events" << endl;
01219   }
01220   //newh->Sumw2();
01221   return newh;
01222 }

vector< Double_t > NuUtilities::RecoBins ( const NuBinningScheme::NuBinningScheme_t  binningScheme  )  [static]

Definition at line 532 of file NuUtilities.cxx.

References GenerateBins(), NuBinningScheme::k0_1GeVTo20, NuBinningScheme::k0_5GeVTo200, NuBinningScheme::kCC0325Std, NuBinningScheme::kCC3flav, NuBinningScheme::kDisplay1GevTo20, NuBinningScheme::kDisplayFD, NuBinningScheme::kDisplayFD2010, NuBinningScheme::kDisplayFD2014, NuBinningScheme::kDisplayFD7e20, NuBinningScheme::kDisplayFDRHC, NuBinningScheme::kDisplayND, Msg::kInfo, NuBinningScheme::kNuMuBar0325Std, NuBinningScheme::kNuMuBar0325Std2, NuBinningScheme::kSterile, NuBinningScheme::kUnknown, and MAXMSG.

Referenced by NuDSTAna::Contamination(), NuDSTAna::ContaminationNQ(), ToFPlotter::CreateCCHistograms(), FoverNHistos::CreateHistograms(), ToFPlotter::CreateNCHistograms(), ToFPlotter::CreateRAFHistograms(), NuHistos::FillMatrixMethodHistos(), NuHistos::FillMatrixMethodNCHistos(), NuDSTAna::MakeResolutionBins(), NuDSTAna::NDOsc(), NuFCExperimentFactory::NuFCExperimentFactory(), NuFCExperimentFactoryNSI::NuFCExperimentFactoryNSI(), NuMatrixMethod::NuMatrixMethod(), NuTransition::NuTransition(), RebinHistogram(), and NuShiftableUnbinnedSpectrum::Spectrum().

00533 {
00534   if (NuBinningScheme::kUnknown == binningScheme ||
00535       NuBinningScheme::k0_5GeVTo200 == binningScheme){
00536     MAXMSG("NuUtilities",Msg::kInfo,5)
00537       << "Binning reco energy in 0.5 GeV up to 200 GeV" << endl;
00538     vector<Double_t> bins;
00539     Double_t lowEdge = 0.0;
00540     Double_t highEdge = 200.0;
00541     Double_t binWidth = 0.5;
00542     for (Double_t binEdge = lowEdge;
00543          binEdge <= highEdge + (binWidth/2.0);
00544          binEdge += binWidth){
00545       bins.push_back(binEdge);
00546     }
00547     return bins;
00548   } else if (NuBinningScheme::kCC0325Std == binningScheme) {
00549     MAXMSG("NuUtilities",Msg::kInfo,5)
00550       << "Reco energy: one bin from 0.0--0.5, "
00551       << "then 0.25 GeV bins up to 200.0 GeV"
00552       << endl;
00553     vector<Double_t> bins;
00554     bins.push_back(0.0);
00555     Double_t lowEdge = 0.5;
00556     Double_t highEdge = 200.0;
00557     Double_t binWidth = 0.25;
00558     for (Double_t binEdge = lowEdge;
00559          binEdge <= highEdge + (binWidth/2.0);
00560          binEdge += binWidth){
00561       bins.push_back(binEdge);
00562     }
00563     return bins;
00564   } else if (NuBinningScheme::kNuMuBar0325Std == binningScheme) {
00565     MAXMSG("NuUtilities",Msg::kInfo,5)
00566       << "Reco energy: 0.5 GeV up to 30 GeV; 2 GeV up to 50 GeV; "
00567       << "1 bin up to 200 GeV"
00568       << endl;
00569     vector<Double_t> bins;
00570     Double_t lowEdge = 0.0;
00571     Double_t highEdge = 30.0;
00572     Double_t binWidth = 0.5;
00573     for (Double_t binEdge = lowEdge;
00574          binEdge < highEdge - (binWidth/2.0);
00575          binEdge += binWidth){
00576       bins.push_back(binEdge);
00577     }
00578     lowEdge = 30.0;
00579     highEdge = 50.0;
00580     binWidth = 2.0;
00581     for (Double_t binEdge = lowEdge;
00582          binEdge < highEdge - (binWidth/2.0);
00583          binEdge += binWidth){
00584       bins.push_back(binEdge);
00585     }
00586     bins.push_back(50.0);
00587     bins.push_back(200.0);
00588     return bins;
00589   } else if (NuBinningScheme::kNuMuBar0325Std2 == binningScheme) {
00590     MAXMSG("NuUtilities",Msg::kInfo,5)
00591       << "Reco Energy: One 0.5 GeV, then 0.25 GeV up to 20 GeV; "
00592       << "1 Gev to 30 GeV,\n 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
00593       << endl;
00594     vector<Double_t> bins;
00595     bins.push_back(0.0);
00596     Double_t lowEdge = 0.5;
00597     Double_t highEdge = 20.0;
00598     Double_t binWidth = 0.25;
00599     for (Double_t binEdge = lowEdge;
00600           binEdge < highEdge -(binWidth/2.0);
00601           binEdge += binWidth) {
00602       bins.push_back(binEdge);
00603     }
00604     lowEdge = 20.0;
00605     highEdge = 30.0;
00606     binWidth = 1.0;
00607     for (Double_t binEdge = lowEdge;
00608           binEdge < highEdge -(binWidth/2.0);
00609           binEdge += binWidth) {
00610       bins.push_back(binEdge);
00611     }
00612     lowEdge = 30.0;
00613     highEdge = 50.0;
00614     binWidth = 2.0;
00615     for (Double_t binEdge = lowEdge;
00616           binEdge < highEdge -(binWidth/2.0);
00617           binEdge += binWidth) {
00618       bins.push_back(binEdge);
00619     }
00620     bins.push_back(50.0);
00621     bins.push_back(200.0);
00622     return bins;
00623   } else if (NuBinningScheme::kCC3flav == binningScheme) {
00624     MAXMSG("NuUtilities",Msg::kInfo,5)
00625       << "Reco Energy: One 1.0 GeV, then 0.25 GeV up to 20 GeV; "
00626       << "1 Gev to 30 GeV,\n 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
00627       << endl;
00628     vector<Double_t> bins;
00629     bins.push_back(0.0);
00630     Double_t lowEdge = 1.0;
00631     Double_t highEdge = 20.0;
00632     Double_t binWidth = 0.25;
00633     for (Double_t binEdge = lowEdge;
00634           binEdge < highEdge -(binWidth/2.0);
00635           binEdge += binWidth) {
00636       bins.push_back(binEdge);
00637     }
00638     lowEdge = 20.0;
00639     highEdge = 30.0;
00640     binWidth = 1.0;
00641     for (Double_t binEdge = lowEdge;
00642           binEdge < highEdge -(binWidth/2.0);
00643           binEdge += binWidth) {
00644       bins.push_back(binEdge);
00645     }
00646     lowEdge = 30.0;
00647     highEdge = 50.0;
00648     binWidth = 2.0;
00649     for (Double_t binEdge = lowEdge;
00650           binEdge < highEdge -(binWidth/2.0);
00651           binEdge += binWidth) {
00652       bins.push_back(binEdge);
00653     }
00654     bins.push_back(50.0);
00655     bins.push_back(200.0);
00656     return bins;
00657   }else if (NuBinningScheme::kSterile == binningScheme) { //A.Devan 8/29/2012
00658     MAXMSG("NuUtilities",Msg::kInfo,5)
00659       << "Reco Energy: 0.05GeV bin between 0 and 5 GeV, then 0.25 GeV up to 20 GeV; "
00660       << "1 Gev to 30 GeV,\n 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
00661       << endl;
00662     vector<Double_t> bins;
00663     bins.push_back(0.0);
00664         Double_t lowEdge = 0.05;
00665         Double_t highEdge = 5.0;
00666         Double_t binWidth = 0.05;
00667         for (Double_t binEdge = lowEdge;
00668                   binEdge < highEdge - (binWidth/2.0);
00669                   binEdge += binWidth) {
00670           bins.push_back(binEdge);
00671         }
00672         lowEdge = 5.0;
00673     highEdge = 20.0;
00674     binWidth = 0.25;
00675     for (Double_t binEdge = lowEdge;
00676           binEdge < highEdge -(binWidth/2.0);
00677           binEdge += binWidth) {
00678       bins.push_back(binEdge);
00679     }
00680     lowEdge = 20.0;
00681     highEdge = 30.0;
00682     binWidth = 1.0;
00683     for (Double_t binEdge = lowEdge;
00684           binEdge < highEdge -(binWidth/2.0);
00685           binEdge += binWidth) {
00686       bins.push_back(binEdge);
00687     }
00688     lowEdge = 30.0;
00689     highEdge = 50.0;
00690     binWidth = 2.0;
00691     for (Double_t binEdge = lowEdge;
00692           binEdge < highEdge -(binWidth/2.0);
00693           binEdge += binWidth) {
00694       bins.push_back(binEdge);
00695     }
00696     bins.push_back(50.0);
00697     bins.push_back(200.0);
00698     return bins;
00699   } else if (NuBinningScheme::kDisplayFD == binningScheme) {
00700     MAXMSG("NuUtilities",Msg::kInfo,5)
00701     << "Reco Energy: 4 GeV Bins up to 20 GeV, then 20, 30, 50, 200"
00702     << endl;
00703     vector<Double_t> bins;
00704     Double_t lowEdge = 0.0;
00705     Double_t highEdge = 20.0;
00706     Double_t binWidth = 4;
00707     for (Double_t binEdge = lowEdge;
00708          binEdge < highEdge -(binWidth/2.0);
00709          binEdge += binWidth) {
00710       bins.push_back(binEdge);
00711     }
00712     bins.push_back(20.0);
00713     bins.push_back(30.0);
00714     bins.push_back(50.0);
00715     bins.push_back(200.0);
00716     return bins;
00717   } else if (NuBinningScheme::kDisplayND == binningScheme) {
00718     MAXMSG("NuUtilities",Msg::kInfo,5)
00719     << "Reco Energy: 1 Gev to 30 GeV, 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
00720     << endl;
00721     vector<Double_t> bins;
00722     Double_t lowEdge = 0.0;
00723     Double_t highEdge = 20.0;
00724     Double_t binWidth = 1.0;
00725     for (Double_t binEdge = lowEdge;
00726          binEdge < highEdge -(binWidth/2.0);
00727          binEdge += binWidth) {
00728       bins.push_back(binEdge);
00729     }
00730     lowEdge = 20.0;
00731     highEdge = 30.0;
00732     binWidth = 2.0;
00733     for (Double_t binEdge = lowEdge;
00734          binEdge < highEdge -(binWidth/2.0);
00735          binEdge += binWidth) {
00736       bins.push_back(binEdge);
00737     }
00738     lowEdge = 30.0;
00739     highEdge = 50.0;
00740     binWidth = 4.0;
00741     for (Double_t binEdge = lowEdge;
00742          binEdge < highEdge -(binWidth/2.0);
00743          binEdge += binWidth) {
00744       bins.push_back(binEdge);
00745     }
00746     bins.push_back(50.0);
00747     bins.push_back(200.0);
00748     return bins;
00749   } else if (NuBinningScheme::kDisplayFDRHC == binningScheme) {
00750       double edges[] = {0, 10, 20, 30,  50, 200};
00751       double steps[] = {  1, 5, 10, 20, 150};
00752       int ng = 5;
00753 
00754       std::vector<Double_t> bins = NuUtilities::GenerateBins(ng, edges, steps);
00755       return bins;
00756   } else if(NuBinningScheme::kDisplayFD2010 == binningScheme){
00757     double edges[] = {0, 1, 5, 10, 20, 30, 50, 200};
00758     double steps[] = {1/*.5*/, .25, .5, 1, 10, 20, 130};
00759     int ng = 7;
00760 
00761     vector<double> bins = NuUtilities::GenerateBins(ng, edges, steps);
00762     return bins;
00763   } else if(NuBinningScheme::kDisplayFD2014 == binningScheme){
00764     double edges[] = {0, 1.5, 10, 20, 30, 50, 200};
00765     double steps[] = {1.5, 0.5, 2, 10, 20, 130};
00766     int ng = 6;
00767 
00768     vector<double> bins = NuUtilities::GenerateBins(ng, edges, steps);
00769     return bins;
00770   } else if(NuBinningScheme::kDisplayFD2014 == binningScheme){
00771     double edges[] = {0, 1.5, 10, 20, 30, 50, 200};
00772     double steps[] = {1.5, 0.5, 2, 10, 20, 130};
00773     int ng = 6;
00774 
00775     vector<double> bins = NuUtilities::GenerateBins(ng, edges, steps);
00776     return bins;
00777   } else if(NuBinningScheme::kDisplayFD7e20 == binningScheme){
00778       double edges[] = {0, 20, 30, 50, 200};
00779       double steps[] = {  2, 5, 20, 150};
00780       int ng = 4;
00781 
00782     vector<double> bins = NuUtilities::GenerateBins(ng, edges, steps);
00783     return bins;
00784   } else if(NuBinningScheme::kDisplay1GevTo20 == binningScheme) {
00785     MAXMSG("NuUtilities", Msg::kInfo, 5)
00786     << "Reco Energy: 0 Gev to 20 GeV in 1 GeV bins, plus a single 20--200 GeV bin"
00787     << endl;
00788     vector<Double_t> bins;
00789     Double_t lowEdge = 0.0;
00790     Double_t highEdge = 20.0;
00791     Double_t binWidth = 1.0;
00792     for (Double_t binEdge = lowEdge;
00793          binEdge < highEdge - (binWidth / 2.0);
00794          binEdge += binWidth) {
00795       bins.push_back(binEdge);
00796     }
00797     bins.push_back(20.0);
00798     bins.push_back(50.0);
00799     bins.push_back(200.0);
00800 
00801     return bins;
00802   } else if (NuBinningScheme::k0_1GeVTo20 == binningScheme) {
00803     MAXMSG("NuUtilities",Msg::kInfo,5)
00804       << "Reco Energy: 0.10 GeV up to 20 GeV; "
00805       << "1 Gev to 30 GeV, 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
00806       << endl;
00807     vector<Double_t> bins;
00808     Double_t lowEdge = 0.0;
00809     Double_t highEdge = 20.0;
00810     Double_t binWidth = 0.1;
00811     for (Double_t binEdge = lowEdge;
00812           binEdge < highEdge -(binWidth/2.0);
00813           binEdge += binWidth) {
00814       bins.push_back(binEdge);
00815     }
00816     lowEdge = 20.0;
00817     highEdge = 30.0;
00818     binWidth = 1.0;
00819     for (Double_t binEdge = lowEdge;
00820           binEdge < highEdge -(binWidth/2.0);
00821           binEdge += binWidth) {
00822       bins.push_back(binEdge);
00823     }
00824     lowEdge = 30.0;
00825     highEdge = 50.0;
00826     binWidth = 2.0;
00827     for (Double_t binEdge = lowEdge;
00828           binEdge < highEdge -(binWidth/2.0);
00829           binEdge += binWidth) {
00830       bins.push_back(binEdge);
00831     }
00832     bins.push_back(50.0);
00833     bins.push_back(200.0);
00834     return bins;
00835   } else {
00836     MAXMSG("NuUtilities",Msg::kInfo,100)
00837     << "!!Invalid binning scheme (" << binningScheme << ") requested.!!" << endl
00838     << "Binning reco energy in 0.5 GeV up to 200 GeV" << endl;
00839     vector<Double_t> bins;
00840     Double_t lowEdge = 0.0;
00841     Double_t highEdge = 200.0;
00842     Double_t binWidth = 0.5;
00843     for (Double_t binEdge = lowEdge;
00844          binEdge <= highEdge + (binWidth/2.0);
00845          binEdge += binWidth){
00846       bins.push_back(binEdge);
00847     }
00848     return bins;
00849   }
00850 } // RecoBins

BeamType::BeamType_t NuUtilities::ReturnConventionsBeamType ( const NuConfig config  )  [static]

Definition at line 342 of file NuUtilities.cxx.

References NuConfig::beamType, ConvertToConventionsBeamType(), NuConfig::detector, NuConfig::intensity, SimFlag::kData, Detector::kFar, Msg::kInfo, MAXMSG, and NuConfig::simFlag.

Referenced by NuBeam::IsGoodSpillAndFillPot().

00343 {
00344   if (Detector::kFar==config.detector && SimFlag::kData==config.simFlag){
00345     MAXMSG("NuUtilities",Msg::kInfo,1)
00346       << "Detected FD data, so not combining intensity with beamType"
00347       << endl;
00348     BeamType::BeamType_t originalBeamType =
00349       static_cast<BeamType::BeamType_t>(config.beamType);
00350     return originalBeamType;
00351   }
00352   else{
00353     return NuUtilities::ConvertToConventionsBeamType(config.beamType,
00354                                                      config.intensity);
00355   }
00356 }

vector< Double_t > NuUtilities::TrueBins ( const NuBinningScheme::NuBinningScheme_t  binningScheme  )  [static]

Definition at line 854 of file NuUtilities.cxx.

References NuBinningScheme::k0_1GeVTo20, NuBinningScheme::k0_5GeVTo200, NuBinningScheme::kCC0325Std, NuBinningScheme::kCC3flav, NuBinningScheme::kDisplayFD, NuBinningScheme::kDisplayND, Msg::kInfo, NuBinningScheme::kNuMuBar0325Std, NuBinningScheme::kNuMuBar0325Std2, NuBinningScheme::kUnknown, and MAXMSG.

Referenced by FoverNHistos::CreateHistograms(), NuHistos::FillMatrixMethodHistos(), NuHistos::FillMatrixMethodNCHistos(), NuPlots::FillTrueFidEnergySpect(), NuFCExperimentFactory::NuFCExperimentFactory(), NuFCExperimentFactoryNSI::NuFCExperimentFactoryNSI(), NuFluxHelper::NuFluxHelper(), NuMatrixMethod::NuMatrixMethod(), and NuTransition::NuTransition().

00855 {
00856   if (NuBinningScheme::kUnknown == binningScheme ||
00857       NuBinningScheme::k0_5GeVTo200 == binningScheme){
00858     MAXMSG("NuUtilities",Msg::kInfo,5)
00859       << "Binning true energy in 0.5 GeV up to 200 GeV" << endl;
00860     vector<Double_t> bins;
00861     Double_t lowEdge = 0.0;
00862     Double_t highEdge = 200.0;
00863     Double_t binWidth = 0.5;
00864     for (Double_t binEdge = lowEdge;
00865          binEdge <= highEdge + (binWidth/2.0);
00866          binEdge += binWidth){
00867       bins.push_back(binEdge);
00868     }
00869     return bins;
00870   }
00871   else if (NuBinningScheme::kCC0325Std == binningScheme){
00872     MAXMSG("NuUtilities",Msg::kInfo,5)
00873       << "True energy: one bin from 0.0--0.5, "
00874       << "then 0.25 GeV bins up to 200.0 GeV"
00875       << endl;
00876     vector<Double_t> bins;
00877     bins.push_back(0.0);
00878     Double_t lowEdge = 0.5;
00879     Double_t highEdge = 200.0;
00880     Double_t binWidth = 0.25;
00881     for (Double_t binEdge = lowEdge;
00882          binEdge <= highEdge + (binWidth/2.0);
00883          binEdge += binWidth){
00884       bins.push_back(binEdge);
00885     }
00886     return bins;
00887   }
00888   else if (NuBinningScheme::kCC3flav == binningScheme) {
00889     MAXMSG("NuUtilities",Msg::kInfo,5)
00890       << "True Energy: 0.05 GeV up to 5 GeV, then 0.25 GeV up to 20 GeV; "
00891       << "1 Gev to 30 GeV,\n 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
00892       << endl;
00893     vector<Double_t> bins;
00894     Double_t lowEdge = 0;
00895     Double_t highEdge = 0.5;
00896     Double_t binWidth = 0.05;
00897     for (Double_t binEdge = lowEdge;
00898           binEdge < highEdge -(binWidth/2.0);
00899           binEdge += binWidth) {
00900       bins.push_back(binEdge);
00901     }
00902     lowEdge = 0.5;
00903     highEdge = 5.0;
00904     binWidth = 0.05;
00905     for (Double_t binEdge = lowEdge;
00906           binEdge < highEdge -(binWidth/2.0);
00907           binEdge += binWidth) {
00908       bins.push_back(binEdge);
00909     }
00910     lowEdge = 5;
00911     highEdge = 20.0;
00912     binWidth = 0.25;
00913     for (Double_t binEdge = lowEdge;
00914           binEdge < highEdge -(binWidth/2.0);
00915           binEdge += binWidth) {
00916       bins.push_back(binEdge);
00917     }
00918     lowEdge = 20.0;
00919     highEdge = 30.0;
00920     binWidth = 1.0;
00921     for (Double_t binEdge = lowEdge;
00922           binEdge < highEdge -(binWidth/2.0);
00923           binEdge += binWidth) {
00924       bins.push_back(binEdge);
00925     }
00926     lowEdge = 30.0;
00927     highEdge = 50.0;
00928     binWidth = 2.0;
00929     for (Double_t binEdge = lowEdge;
00930           binEdge < highEdge -(binWidth/2.0);
00931           binEdge += binWidth) {
00932       bins.push_back(binEdge);
00933     }
00934     bins.push_back(50.0);
00935     bins.push_back(200.0);
00936     return bins;
00937   }
00938   else if (NuBinningScheme::kDisplayFD == binningScheme) {
00939     MAXMSG("NuUtilities",Msg::kInfo,5)
00940     << "True Energy: 4 GeV Bins up to 20 GeV, then 20, 30, 50, 200"
00941     << endl;
00942     vector<Double_t> bins;
00943     Double_t lowEdge = 0.0;
00944     Double_t highEdge = 20.0;
00945     Double_t binWidth = 4;
00946     for (Double_t binEdge = lowEdge;
00947          binEdge < highEdge -(binWidth/2.0);
00948          binEdge += binWidth) {
00949       bins.push_back(binEdge);
00950     }
00951     bins.push_back(20.0);
00952     bins.push_back(30.0);
00953     bins.push_back(50.0);
00954     bins.push_back(200.0);
00955     return bins;
00956   }
00957   else if (NuBinningScheme::kNuMuBar0325Std == binningScheme){
00958     MAXMSG("NuUtilities",Msg::kInfo,5)
00959       << "True energy: 0.5 GeV up to 30 GeV; 2 GeV up to 50 GeV; "
00960       << "1 bin up to 200 GeV"
00961       << endl;
00962     vector<Double_t> bins;
00963     Double_t lowEdge = 0.0;
00964     Double_t highEdge = 30.0;
00965     Double_t binWidth = 0.5;
00966     for (Double_t binEdge = lowEdge;
00967          binEdge < highEdge - (binWidth/2.0);
00968          binEdge += binWidth){
00969       bins.push_back(binEdge);
00970     }
00971     lowEdge = 30.0;
00972     highEdge = 50.0;
00973     binWidth = 2.0;
00974     for (Double_t binEdge = lowEdge;
00975          binEdge < highEdge - (binWidth/2.0);
00976          binEdge += binWidth){
00977       bins.push_back(binEdge);
00978     }
00979     bins.push_back(50.0);
00980     bins.push_back(200.0);
00981     return bins;
00982   }
00983   else if (NuBinningScheme::kNuMuBar0325Std2 == binningScheme) {
00984     MAXMSG("NuUtilities",Msg::kInfo,5)
00985       << "True Energy: One 0.5 GeV, then 0.25 GeV up to 20 GeV; "
00986       << "1 Gev to 30 GeV,\n 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
00987       << endl;
00988     vector<Double_t> bins;
00989     bins.push_back(0.0);
00990     Double_t lowEdge = 0.5;
00991     Double_t highEdge = 20.0;
00992     Double_t binWidth = 0.25;
00993     for (Double_t binEdge = lowEdge;
00994           binEdge < highEdge -(binWidth/2.0);
00995           binEdge += binWidth) {
00996       bins.push_back(binEdge);
00997     }
00998     lowEdge = 20.0;
00999     highEdge = 30.0;
01000     binWidth = 1.0;
01001     for (Double_t binEdge = lowEdge;
01002           binEdge < highEdge -(binWidth/2.0);
01003           binEdge += binWidth) {
01004       bins.push_back(binEdge);
01005     }
01006     lowEdge = 30.0;
01007     highEdge = 50.0;
01008     binWidth = 2.0;
01009     for (Double_t binEdge = lowEdge;
01010           binEdge < highEdge -(binWidth/2.0);
01011           binEdge += binWidth) {
01012       bins.push_back(binEdge);
01013     }
01014     bins.push_back(50.0);
01015     bins.push_back(200.0);
01016     return bins;
01017   }
01018   else if (NuBinningScheme::kDisplayFD == binningScheme) {
01019     MAXMSG("NuUtilities",Msg::kInfo,5)
01020     << "True Energy: 4 GeV Bins up to 20 GeV, then 20, 30, 50, 200"
01021     << endl;
01022     vector<Double_t> bins;
01023     Double_t lowEdge = 0.0;
01024     Double_t highEdge = 20.0;
01025     Double_t binWidth = 4;
01026     for (Double_t binEdge = lowEdge;
01027          binEdge < highEdge -(binWidth/2.0);
01028          binEdge += binWidth) {
01029       bins.push_back(binEdge);
01030     }
01031     bins.push_back(20.0);
01032     bins.push_back(30.0);
01033     bins.push_back(50.0);
01034     bins.push_back(200.0);
01035     return bins;
01036   }
01037   else if (NuBinningScheme::kDisplayND == binningScheme) {
01038     MAXMSG("NuUtilities",Msg::kInfo,5)
01039     << "Reco Energy: 1 Gev to 30 GeV, 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
01040     << endl;
01041     vector<Double_t> bins;
01042     Double_t lowEdge = 0.0;
01043     Double_t highEdge = 20.0;
01044     Double_t binWidth = 1.0;
01045     for (Double_t binEdge = lowEdge;
01046          binEdge < highEdge -(binWidth/2.0);
01047          binEdge += binWidth) {
01048       bins.push_back(binEdge);
01049     }
01050     lowEdge = 20.0;
01051     highEdge = 30.0;
01052     binWidth = 2.0;
01053     for (Double_t binEdge = lowEdge;
01054          binEdge < highEdge -(binWidth/2.0);
01055          binEdge += binWidth) {
01056       bins.push_back(binEdge);
01057     }
01058     lowEdge = 30.0;
01059     highEdge = 50.0;
01060     binWidth = 4.0;
01061     for (Double_t binEdge = lowEdge;
01062          binEdge < highEdge -(binWidth/2.0);
01063          binEdge += binWidth) {
01064       bins.push_back(binEdge);
01065     }
01066     bins.push_back(50.0);
01067     bins.push_back(200.0);
01068     return bins;
01069   } else if (NuBinningScheme::k0_1GeVTo20 == binningScheme) {
01070     MAXMSG("NuUtilities",Msg::kInfo,5)
01071       << "True Energy: 0.10 GeV up to 20 GeV; "
01072       << "1 Gev to 30 GeV, 2 GeV from 30 to 50 GeV; One 150 GeV bin to 200 GeV"
01073       << endl;
01074     vector<Double_t> bins;
01075     Double_t lowEdge = 0.0;
01076     Double_t highEdge = 20.0;
01077     Double_t binWidth = 0.1;
01078     for (Double_t binEdge = lowEdge;
01079           binEdge < highEdge -(binWidth/2.0);
01080           binEdge += binWidth) {
01081       bins.push_back(binEdge);
01082     }
01083     lowEdge = 20.0;
01084     highEdge = 30.0;
01085     binWidth = 1.0;
01086     for (Double_t binEdge = lowEdge;
01087           binEdge < highEdge -(binWidth/2.0);
01088           binEdge += binWidth) {
01089       bins.push_back(binEdge);
01090     }
01091     lowEdge = 30.0;
01092     highEdge = 50.0;
01093     binWidth = 2.0;
01094     for (Double_t binEdge = lowEdge;
01095           binEdge < highEdge -(binWidth/2.0);
01096           binEdge += binWidth) {
01097       bins.push_back(binEdge);
01098     }
01099     bins.push_back(50.0);
01100     bins.push_back(200.0);
01101     return bins;
01102  } else {
01103     MAXMSG("NuUtilities",Msg::kInfo,100)
01104     << "!!Invalid binning scheme requested.!!" << endl
01105     << "Binning true energy in 0.5 GeV up to 200 GeV" << endl;
01106     vector<Double_t> bins;
01107     Double_t lowEdge = 0.0;
01108     Double_t highEdge = 200.0;
01109     Double_t binWidth = 0.5;
01110     for (Double_t binEdge = lowEdge;
01111          binEdge <= highEdge + (binWidth/2.0);
01112          binEdge += binWidth){
01113       bins.push_back(binEdge);
01114     }
01115     return bins;
01116   }
01117 }

static Double_t NuUtilities::ValueAt ( TH1 *  h,
Double_t  val 
) [inline, static]

Definition at line 65 of file NuUtilities.h.

Referenced by NuTransSME::Reweight(), and NuTransition::Reweight().

00066     {return h->GetBinContent(h->GetXaxis()->FindFixBin(val));};

void NuUtilities::WriteHistosToFile ( const std::pair< NuMatrixSpectrum, NuMatrixSpectrum > &  vHistos,
TFile &  file 
) [static]

Definition at line 1154 of file NuUtilities.cxx.

01156 {
01157   file.cd();
01158   vHistos.first.Spectrum()->Write();
01159   vHistos.second.Spectrum()->Write();
01160 }

void NuUtilities::WriteHistosToFile ( const std::vector< NuMatrixSpectrum > &  vHistos,
TFile &  file,
bool  Quiet = false 
) [static]

Definition at line 1137 of file NuUtilities.cxx.

01139 {
01140   file.cd();
01141   for (vector<NuMatrixSpectrum>::const_iterator it = vHistos.begin();
01142        it != vHistos.end();
01143        ++it){
01144     (*it).Spectrum()->Write();
01145     if(!Quiet) cout << "Written " << (*it).GetName() << endl;
01146   }
01147 
01148   cout << "Written " << vHistos.size() << " histograms to " 
01149        << file.GetName() << endl;
01150 
01151 }

void NuUtilities::WriteHistosToFile ( const std::vector< TH1D > &  vHistos,
TFile &  file,
bool  Quiet = false 
) [static]

Definition at line 1120 of file NuUtilities.cxx.

01122 {
01123   file.cd();
01124   for (vector<TH1D>::const_iterator it = vHistos.begin();
01125        it != vHistos.end();
01126        ++it){
01127     (*it).Write();
01128     if(!Quiet) cout << "Written " << (*it).GetName() << endl;
01129   }
01130 
01131   cout << "Written " << vHistos.size() << " histograms to " 
01132        << file.GetName() << endl;
01133 
01134 }


Member Data Documentation

VldTimeStamp NuUtilities::vtsRunIEnd [static]

Definition at line 101 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunIIEnd [static]

Definition at line 102 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunIIIEnd [static]

Definition at line 103 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunIVEnd [static]

Definition at line 104 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunIXEnd [static]

Definition at line 109 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunVEnd [static]

Definition at line 105 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunVIEnd [static]

Definition at line 106 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunVIIEnd [static]

Definition at line 107 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunVIIIEnd [static]

Definition at line 108 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().

VldTimeStamp NuUtilities::vtsRunXEnd [static]

Definition at line 110 of file NuUtilities.h.

Referenced by NuExtraction::ExtractBeamInfoDB(), and NuAnalysis::ExtractConfig().


The documentation for this class was generated from the following files:
Generated on Mon Mar 2 00:38:00 2015 for loon by  doxygen 1.4.7