Zbeam Class Reference

#include <Zbeam.h>

List of all members.

Classes

struct  ZbeamData

Public Types

typedef struct Zbeam::ZbeamData ZbeamData_t

Public Member Functions

 Zbeam (const char *dir="")
virtual ~Zbeam ()
void SetReweightConfig (std::string cnf)
double GetWeight (ZbeamData_t zbmdata, BeamSys::BeamSys_t sys, double nSigma, bool is120GeVband=false)
double GetWeight (ZbeamData_t zbmdata, std::map< BeamSys::BeamSys_t, double > nSigma)
double GetFarOverNearError (ZbeamData_t zbmdata, BeamSys::BeamSys_t sys, double nSigma)
double GetWeight (int iDet, int iBeam, int nEffect, double nSigma, int ntype, double Enu)
double GetWeight (int iDet, int iBeam, std::vector< double > nSigma, int ntype, double Enu)
double GetWeight (int iDet, int iBeam, int nEffect, double nSigma, double Enu)
void MakeHistogram (std::string name, std::vector< double > vec)
void SaveHistogram (std::string name, std::string config, std::string file)
void DeleteHistogram (std::string name, std::string config, std::string file)
void FillVector (TH1D *hist, int ntype, int eff, int beam, int det)
TFile * GetBeamSysFile ()

Private Types

typedef std::map< int, double * > EffMap_t
typedef std::map< int, EffMap_tBeamMap_t
typedef std::map< int, BeamMap_tDetMap_t

Private Attributes

TFile * fBeamSysFile
std::string fCurrentConfig
std::string fTopDir
std::map< int, DetMap_tfBeamSysMap
std::map< std::string, TH1D * > fHistMap

Detailed Description

Definition at line 152 of file Zbeam.h.


Member Typedef Documentation

typedef std::map<int , EffMap_t > Zbeam::BeamMap_t [private]

Definition at line 207 of file Zbeam.h.

typedef std::map<int , BeamMap_t > Zbeam::DetMap_t [private]

Definition at line 208 of file Zbeam.h.

typedef std::map<int , double* > Zbeam::EffMap_t [private]

Definition at line 206 of file Zbeam.h.


Constructor & Destructor Documentation

Zbeam::Zbeam ( const char *  dir = ""  ) 

Definition at line 149 of file Zbeam.cxx.

References base, gSystem(), Msg::kDebug, Msg::kFatal, and MSG.

00150 {
00151   // fTopDir is where the beamsys_hists.root is
00152   fTopDir=dir; // user may set location of input data
00153   if(fTopDir=="") { // by default, this code looks in a standard place
00154     fTopDir="MCReweight/data";
00155     std::string base="";
00156     base=getenv("SRT_PRIVATE_CONTEXT");
00157     if (base!="" && base!=".")
00158       {
00159         // check if directory exists in SRT_PRIVATE_CONTEXT
00160         std::string path = base + "/" + fTopDir;
00161         void *dir_ptr = gSystem -> OpenDirectory(path.c_str());
00162         if(!dir_ptr) base=getenv("SRT_PUBLIC_CONTEXT"); // if it doesn't exist then use SRT_PUBLIC_CONTEXT
00163       }
00164     else base=getenv("SRT_PUBLIC_CONTEXT");
00165     
00166     if(base=="") {
00167       MSG("MCReweight",Msg::kFatal)<<"No SRT_PUBLIC_CONTEXT set"<<std::endl;
00168       assert(false);
00169     }
00170     fTopDir = base+ "/" + fTopDir;
00171   }
00172     MSG("Zbeam",Msg::kDebug) <<"Zbeam getting hists from: "<<fTopDir<<"/beamsys_hists.root"<<std::endl;
00173   
00174   SetReweightConfig("DetXs"); //default config
00175 }

Zbeam::~Zbeam (  )  [virtual]

Definition at line 177 of file Zbeam.cxx.

00177               {
00178 
00179 }


Member Function Documentation

void Zbeam::DeleteHistogram ( std::string  name,
std::string  config,
std::string  file 
)

Definition at line 440 of file Zbeam.cxx.

00441 {
00442   TFile* f=new TFile(filename.c_str(),"UPDATE");
00443   f->cd(config.c_str());
00444   f->rmdir(name.c_str());
00445   f->Close();
00446   delete f;
00447   f=0;
00448 }

void Zbeam::FillVector ( TH1D *  hist,
int  ntype,
int  eff,
int  beam,
int  det 
)

Definition at line 450 of file Zbeam.cxx.

References fBeamSysMap, and find().

Referenced by SetReweightConfig().

00451 {
00452   EffMap_t map1;
00453   BeamMap_t map2;
00454   DetMap_t map3;
00455 
00456   const int NBINS=hist->GetNbinsX();
00457   double* vec=new double[NBINS];
00458 
00459   for (int i=1; i<hist->GetNbinsX()+1;i++) vec[i-1]=hist->GetBinContent(i);
00460 
00461   map1.insert(EffMap_t::value_type(eff, vec));
00462   map2.insert(BeamMap_t::value_type(beam, map1));
00463   map3.insert(DetMap_t::value_type(det, map2));
00464 
00465   std::map<int, DetMap_t>::iterator bsmit = fBeamSysMap.find(ntype);
00466 
00467   if (bsmit == fBeamSysMap.end()) {
00468     fBeamSysMap.insert(std::map<int, DetMap_t>::value_type(ntype, map3));
00469   } else {
00470     DetMap_t::iterator dmit = (bsmit->second).find(det);
00471     if (dmit == (bsmit->second).end() ) {
00472       ((*bsmit).second).insert(DetMap_t::value_type(det,map2));
00473     } else {
00474       BeamMap_t::iterator bmit = (dmit->second).find(beam);
00475       if (bmit == (dmit->second).end() ) {
00476         ((*dmit).second).insert(BeamMap_t::value_type(beam,map1));
00477       } else {
00478         EffMap_t::iterator effit = (bmit->second).find(eff);
00479         if (effit == (bmit->second).end() ) {
00480           ((*bmit).second).insert(EffMap_t::value_type(eff,vec));
00481         } else {
00482           cout << "Already loaded vector for: " 
00483                << "beam "<<beam<< " det "<<det<<" nu "<<ntype<<" eff "<<eff<<endl;
00484           return;
00485         }
00486       }
00487     }
00488   } 
00489 }

TFile* Zbeam::GetBeamSysFile (  )  [inline]

Definition at line 199 of file Zbeam.h.

References fBeamSysFile.

00199 {return fBeamSysFile;}

double Zbeam::GetFarOverNearError ( ZbeamData_t  zbmdata,
BeamSys::BeamSys_t  sys,
double  nSigma 
) [inline]

Definition at line 170 of file Zbeam.h.

References Zbeam::ZbeamData::detector, GetWeight(), and Detector::kUnknown.

00171   {
00172     zbmdata.detector = Detector::kUnknown;
00173     return GetWeight(zbmdata, sys, nSigma);
00174   };

double Zbeam::GetWeight ( int  iDet,
int  iBeam,
int  nEffect,
double  nSigma,
double  Enu 
) [inline]

Definition at line 186 of file Zbeam.h.

References GetWeight().

Referenced by GetWeight().

00187   { return GetWeight(iDet, iBeam, nEffect, nSigma, 56, Enu);} 

double Zbeam::GetWeight ( int  iDet,
int  iBeam,
std::vector< double >  nSigma,
int  ntype,
double  Enu 
)

Definition at line 360 of file Zbeam.cxx.

References Zbeam::ZbeamData::beam, Zbeam::ZbeamData::detector, BeamType::FromZarko(), GetWeight(), Detector::kUnknown, Zbeam::ZbeamData::ntype, and Zbeam::ZbeamData::true_enu.

00362 {
00363   ZbeamData_t zbm;
00364   zbm.beam = BeamType::FromZarko(iBeam);
00365   zbm.detector = (iDet != 3) ?  Detector::EDetector(iDet) : Detector::kUnknown;
00366   zbm.ntype = ntype;
00367   zbm.true_enu = Enu;
00368 
00369   std::map<BeamSys::BeamSys_t, double> sgm;
00370   for (unsigned int nEffect=0; nEffect<nSigma.size();nEffect++)
00371     {
00372       BeamSys::BeamSys_t sys = BeamSys::EBeamSys(nEffect+1);
00373       const std::pair<BeamSys::BeamSys_t, double> sp(sys, nSigma[nEffect]);
00374       sgm.insert(sp);
00375     }
00376   return GetWeight(zbm, sgm);
00377 }

double Zbeam::GetWeight ( int  iDet,
int  iBeam,
int  nEffect,
double  nSigma,
int  ntype,
double  Enu 
)

Definition at line 347 of file Zbeam.cxx.

References Zbeam::ZbeamData::beam, Zbeam::ZbeamData::detector, BeamType::FromZarko(), GetWeight(), Detector::kUnknown, Zbeam::ZbeamData::ntype, and Zbeam::ZbeamData::true_enu.

00349 {
00350   ZbeamData_t zbm;
00351   zbm.beam = BeamType::FromZarko(iBeam);
00352   zbm.detector = (iDet != 3) ?  Detector::EDetector(iDet) : Detector::kUnknown;
00353   zbm.ntype = ntype;
00354   zbm.true_enu = Enu; 
00355   BeamSys::BeamSys_t sys = BeamSys::EBeamSys(nEffect);
00356 
00357   return GetWeight(zbm,sys, nSigma);
00358 }

double Zbeam::GetWeight ( ZbeamData_t  zbmdata,
std::map< BeamSys::BeamSys_t, double >  nSigma 
)

Definition at line 334 of file Zbeam.cxx.

References GetWeight().

00335 { 
00336   double wgh=1.;
00337   std::map<BeamSys::BeamSys_t, double>::const_iterator bsit;
00338 
00339   for (bsit=nSigma.begin();bsit != nSigma.end() ;bsit++)
00340     {
00341       wgh*=GetWeight(zbmdata,bsit->first, bsit->second); 
00342     }
00343   cout<<"zbeam weight= "<<wgh<<endl;
00344   return wgh; 
00345 } 

double Zbeam::GetWeight ( ZbeamData_t  zbmdata,
BeamSys::BeamSys_t  sys,
double  nSigma,
bool  is120GeVband = false 
)

Definition at line 231 of file Zbeam.cxx.

References Zbeam::ZbeamData::beam, det, Zbeam::ZbeamData::detector, fBeamSysMap, find(), BeamSys::kBaffleScraping, BeamSys::kBeamWidth, BeamSys::kEndOfList, BeamSys::kFirst, BeamSys::kHadProdAfter, BeamSys::kHadProdBefore, BeamSys::kHorn1Offset, BeamSys::kHornIDist, BeamSys::kTotalAfter, BeamSys::kTotalBefore, BeamSys::kTotalFocusing, Zbeam::ZbeamData::ntype, and Zbeam::ZbeamData::true_enu.

Referenced by SKZPWeightCalculator::GetBeamWeight(), GetFarOverNearError(), SKZPWeightCalculator::GetFluxError(), SKZPWeightCalculator::GetRunPeriodWeight(), and GetWeight().

00232 {
00233 
00234   int ntype  = zbmdata.ntype;
00235   double Enu = zbmdata.true_enu;
00236   Detector::Detector_t det = zbmdata.detector;
00237   BeamType::BeamType_t beam = zbmdata.beam;
00238 
00239   //   cout<<"beam type = "<<beam<<endl;
00240   //note: the following bin number is not equal to the ROOT histogram bin number. This bin number is being used to witha vector that was filled from a histogram. So here for the vector the first bin = 0 and not 1. See the FillVector method.
00241   int bin;
00242   bin=int(2*Enu);
00243   if (bin<0) bin=0;
00244   if ( (beam==44||beam==47) && !is120GeVband){
00245      if (bin>59) bin=59;
00246      } 
00247      else {
00248        if (bin>99) bin=99;
00249      }
00250      //MINOS+ hists only go to 30GeV, MINOS to 50GeV - ANNA
00251   if ( is120GeVband && (bin>179) ) bin=179; //Setting limit at bin 179, [89.5-90) GeV, because the last 30 GeV has very ropey stats
00252 
00253 
00254 
00255   if (ntype == 14) ntype = 56;
00256   else if (ntype == -14) ntype = 55;
00257   else if (ntype == 12) ntype = 53;
00258   else if (ntype == -12) ntype = 52;
00259 
00260   double weight = 1. ;
00261   double totError=0.;
00262   BeamSys::BeamSys_t first=sys; 
00263   BeamSys::BeamSys_t last=sys;
00264   if (sys == BeamSys::kTotalBefore || sys == BeamSys::kTotalAfter) //Total Error
00265     {
00266       first = BeamSys::kFirst;
00267       last = BeamSys::kEndOfList; 
00268     }
00269   else if(sys == BeamSys::kTotalFocusing)
00270     {
00271       first = BeamSys::kFirst;
00272       last = BeamSys::kHornIDist;
00273     } 
00274 
00275   std::map<int, DetMap_t>::iterator bsmit = fBeamSysMap.find(ntype);   
00276   DetMap_t::iterator dmit;
00277   BeamMap_t::iterator bmit;
00278   bool beam_it_exists = true;
00279   if (bsmit == fBeamSysMap.end()) {
00280     weight = 0.; 
00281     beam_it_exists = false;
00282   } else {
00283     dmit = (bsmit->second).find(det);
00284     if (dmit == (bsmit->second).end() ) {
00285       weight = 0.; 
00286       beam_it_exists = false;
00287     } else {
00288       bmit = (dmit->second).find(beam);
00289       if (bmit == (dmit->second).end() ) {
00290         weight = 0.; 
00291         beam_it_exists = false;
00292       } 
00293     }
00294   }
00295 
00296   if (beam_it_exists) {
00297     for (int eff=first;eff<=last;eff++) {    
00298       if (BeamSys::EBeamSys(eff) != BeamSys::kTotalBefore && 
00299           BeamSys::EBeamSys(eff) != BeamSys::kTotalAfter &&
00300           BeamSys::EBeamSys(eff) != BeamSys::kHadProdBefore && 
00301           BeamSys::EBeamSys(eff) != BeamSys::kHadProdAfter &&
00302           BeamSys::EBeamSys(eff) != BeamSys::kTotalFocusing )
00303         {
00304           EffMap_t::iterator effit = (bmit->second).find(eff);
00305           weight =  ( effit == (bmit->second).end() ) ?  0. : ((*effit).second)[bin];
00306           totError += weight*weight;
00307         }
00308     } // loop over all effects 
00309     
00310     if ((sys == BeamSys::kTotalBefore || sys == BeamSys::kHadProdBefore )) {
00311       EffMap_t::iterator effit = (bmit->second).find(BeamSys::kHadProdBefore);
00312       weight =  ( effit == (bmit->second).end() ) ?  0. : ((*effit).second)[bin];
00313       totError += weight*weight;
00314     }
00315     if ((sys == BeamSys::kTotalAfter || sys == BeamSys::kHadProdAfter )) {
00316       EffMap_t::iterator effit = (bmit->second).find(BeamSys::kHadProdAfter);
00317       weight =  ( effit == (bmit->second).end() ) ?  0. : ((*effit).second)[bin];
00318       totError += weight*weight;
00319     }
00320   }
00321 
00322   totError=sqrt(totError);
00323   
00324   if (!(sys == BeamSys::kTotalBefore || sys == BeamSys::kTotalAfter 
00325         || sys == BeamSys::kTotalFocusing))  totError *= ( (weight>0.) - (weight<0.) );
00326   
00327   double nS = (sys == BeamSys::kBeamWidth) ? nSigma + 1. : nSigma; // by default reweights pME, pHE and lowint 
00328   //some of the effects are symmetric and only abs(nS) makes sense
00329   totError = ((sys == BeamSys::kHorn1Offset) ||
00330               (sys == BeamSys::kBaffleScraping)) ? totError*fabs(nS)+1.0 :  totError*nS+1.; 
00331   return totError;  
00332 } 

void Zbeam::MakeHistogram ( std::string  name,
std::vector< double >  vec 
)

Definition at line 379 of file Zbeam.cxx.

References BeamType::AsString(), BeamSys::AsString(), Detector::AsString(), det, fHistMap, BeamType::kEndOfList, BeamSys::kEndOfList, Detector::kFar, BeamSys::kFirst, Detector::kNear, BeamType::kUnknown, NueConvention::nue, and NueConvention::numu.

00380 {
00381   //Function creates a histogram using a vector vec (assumes 0.5GeV bins)
00382   //Check if the name contains the detector, beam and effect
00383   //This doesn't guarantee that the name is correct. 
00384   //The name should look like this NuMu_Horn1Offset_L010z185i_Near ...
00385   bool validname=false;
00386   std::string numu="NuMu";
00387   std::string numubar="NuMuBar";
00388   std::string nue="NuE";
00389   std::string nuebar="NuEBar";
00390   
00391   for (int det=Detector::kNear;det<=Detector::kFar;det++) {
00392     for (int beam=BeamType::kUnknown;beam<BeamType::kEndOfList;beam++) {
00393       for (int eff=BeamSys::kFirst;eff<BeamSys::kEndOfList;eff++) {
00394         if (((name.find(Detector::AsString(Detector::Detector_t(det))) != string::npos)||name.find("FN") != string::npos) &&
00395             (name.find(BeamType::AsString(BeamType::BeamType_t(beam))) != string::npos) &&
00396             (name.find(BeamSys::AsString(BeamSys::BeamSys_t(eff))) != string::npos) &&
00397             ((name.find(numu) != string::npos) || (name.find(numubar) != string::npos) ||
00398              (name.find(nue) != string::npos) || (name.find(nuebar) != string::npos))) {
00399           validname=true;
00400         } 
00401       }
00402     }
00403   }
00404   if (!validname) {
00405     cout<<"Invalid histogram name! "<<name <<endl;
00406     cout<<"Name should include beam (as in BeamType), detector (as in Detector)" 
00407         <<" and beam systematics type"<<endl;
00408     return;
00409   }
00410   cout<<"Creating "<<name<<" histogram"<<endl;
00411   unsigned int nbins=vec.size()-1;
00412   double maxx=0.5*nbins;
00413 
00414   TH1D* hist=new TH1D(name.c_str(),name.c_str(),nbins,0.,maxx);
00415   for (unsigned int i=1;i<nbins+1;i++)
00416     hist->SetBinContent(i,vec[i-1]);
00417   
00418   fHistMap[name]=hist;
00419 }

void Zbeam::SaveHistogram ( std::string  name,
std::string  config,
std::string  file 
)

Definition at line 421 of file Zbeam.cxx.

References fHistMap, and it.

00422 {
00423   map<std::string, TH1D* >::iterator it = fHistMap.find(name);
00424   if ( it == fHistMap.end() ) {
00425     cout << name <<" does not exist!"<<endl;
00426     return;
00427   }
00428   TFile* f=new TFile(filename.c_str(),"UPDATE");
00429   if (!f->GetListOfKeys()->Contains(config.c_str())) {
00430     cout << "Creating directory for " << config << " configuration" << endl;
00431     f->mkdir(config.c_str());
00432   }
00433   f->cd(config.c_str());
00434   fHistMap[name]->Write();
00435   f->Close();
00436   delete f;
00437   f=0;
00438 }

void Zbeam::SetReweightConfig ( std::string  cnf  ) 

Definition at line 181 of file Zbeam.cxx.

References bfld::AsString(), det, fBeamSysFile, fBeamSysMap, fCurrentConfig, FillVector(), Form(), fTopDir, BeamType::kEndOfList, BeamSys::kEndOfList, Detector::kFar, Msg::kFatal, BeamSys::kFirst, Detector::kNear, Detector::kUnknown, BeamType::kUnknown, and MSG.

Referenced by SKZPWeightCalculator::Config(), SKZPWeightCalculator::GetFluxError(), and NuZBeamReweight::GetSKZPWeightCalculator2018().

00182 {
00183   if ( cnf == fCurrentConfig) return;
00184   TDirectory *save = gDirectory;
00185 
00186    fBeamSysFile = new TFile(Form("%s/beamsys_%s.root",fTopDir.c_str(),cnf.c_str()));
00187 
00188   if (fBeamSysFile->IsZombie()) {
00189     MSG("Zbeam",Msg::kFatal)<<"Don't recognize "
00190                             <<cnf
00191                             <<" configuration."
00192                             <<" Config set to :"
00193                             <<fCurrentConfig<<endl;
00194     
00195     return;
00196   } 
00197   cout<<"The config version is: "<<cnf<<" and current config "<<fCurrentConfig<<endl;
00198   
00199   fBeamSysMap.clear();
00200   fCurrentConfig = cnf;
00201   std::string nus[]={"NuMu","NuMuBar","NuE","NuEBar"};
00202   int ntype[]    ={    56,       55  , 53  , 52};
00203   for (int inu=0;inu<4;inu++) {
00204     for (int eff=BeamSys::kFirst;eff<BeamSys::kEndOfList;eff++) {
00205       for (int beam=BeamType::kUnknown;beam<BeamType::kEndOfList;beam++) {
00206         for (int det=Detector::kNear;det<=Detector::kFar;det++) {
00207           std::string hname=nus[inu]+"_"+
00208             BeamSys::AsString(BeamSys::BeamSys_t(eff))+"_"+
00209             BeamType::AsString(BeamType::BeamType_t(beam))+"_"+
00210             Detector::AsString(Detector::Detector_t(det));
00211           TH1D* hist=dynamic_cast<TH1D*> (fBeamSysFile->Get(hname.c_str()));
00212           //      cout<<eff<<" hname is "<<hname<<endl;
00213           if (hist) FillVector(hist,ntype[inu],eff,beam,det);
00214         }
00215         //for far/near error bar internally use Detector::kUnknown
00216         std::string hname=nus[inu]+"_"+
00217           BeamSys::AsString(BeamSys::BeamSys_t(eff))+"_"+
00218           BeamType::AsString(BeamType::BeamType_t(beam))+"_FN";
00219         TH1D* hist=dynamic_cast<TH1D*> (fBeamSysFile->Get(hname.c_str()));
00220         if (hist) FillVector(hist,ntype[inu],eff,beam,Detector::kUnknown);
00221       }
00222     }
00223   }
00224   
00225   fBeamSysFile->Close();
00226   delete fBeamSysFile;
00227   fBeamSysFile = 0;
00228   save->cd();
00229 }


Member Data Documentation

TFile* Zbeam::fBeamSysFile [private]

Definition at line 202 of file Zbeam.h.

Referenced by GetBeamSysFile(), and SetReweightConfig().

std::map<int, DetMap_t > Zbeam::fBeamSysMap [private]

Definition at line 209 of file Zbeam.h.

Referenced by FillVector(), GetWeight(), and SetReweightConfig().

std::string Zbeam::fCurrentConfig [private]

Definition at line 203 of file Zbeam.h.

Referenced by SetReweightConfig().

std::map<std::string, TH1D* > Zbeam::fHistMap [private]

Definition at line 211 of file Zbeam.h.

Referenced by MakeHistogram(), and SaveHistogram().

std::string Zbeam::fTopDir [private]

Definition at line 204 of file Zbeam.h.

Referenced by SetReweightConfig().


The documentation for this class was generated from the following files:

Generated on 15 Nov 2018 for loon by  doxygen 1.6.1