NuCutImps::CC2014NewCoilCutResolution Class Reference

#include <NuCutImps.h>

Inheritance diagram for NuCutImps::CC2014NewCoilCutResolution:
NuCut

List of all members.

Public Member Functions

 CC2014NewCoilCutResolution (const NuPlots *plots, TString ana, int quantile)
virtual ~CC2014NewCoilCutResolution ()
virtual Bool_t InFidVol (const NuEvent &nu) const
void Preselection (const NuEvent &nu)
void Selection (const NuEvent &nu)

Protected Member Functions

void Resolution (const NuEvent &nu)
void SetCutFile (const char *filename)

Protected Attributes

TTree * fResCuts_r1
TTree * fResCuts_r2
TTree * fResCuts_r3
TTree * fResCuts_phe
int fQuantile

Private Attributes

CC2014NewCoilCut fCC2014NewCoilCut

Detailed Description

Definition at line 314 of file NuCutImps.h.


Constructor & Destructor Documentation

NuCutImps::CC2014NewCoilCutResolution::CC2014NewCoilCutResolution ( const NuPlots plots,
TString  ana,
int  quantile 
)

Definition at line 2061 of file NuCutImps.cxx.

References fQuantile, NuLibrary::general, NuGeneral::GetReleaseDirToUse(), NuLibrary::Instance(), NuCuts::kCC2014NewCoilCut, NuCut::SetAnaVersion(), SetCutFile(), and NuCut::SetFidVol().

02064     : NuCut(ana, plots), fResCuts_r1(0), fResCuts_r2(0),
02065       fResCuts_r3(0), fResCuts_phe(0), fQuantile(quantile)
02066   {
02067     assert(fQuantile >= 0 && fQuantile <= 4);
02068 
02069     SetFidVol("cc2008");
02070     SetAnaVersion(NuCuts::kCC2014NewCoilCut);
02071 
02072     const NuGeneral& gen = NuLibrary::Instance().general;
02073     const TString releaseDir = gen.GetReleaseDirToUse("NtupleUtils");
02074     const TString fullPath = releaseDir+"/NtupleUtils/data/cuts_numu.root";
02075     SetCutFile(fullPath);
02076   }

NuCutImps::CC2014NewCoilCutResolution::~CC2014NewCoilCutResolution (  )  [virtual]

Definition at line 2078 of file NuCutImps.cxx.

02079   {
02080   }


Member Function Documentation

Bool_t NuCutImps::CC2014NewCoilCutResolution::InFidVol ( const NuEvent nu  )  const [virtual]

Fiducial volume calculation function. This is the basic 'infid' derived version, and should be overridden/passed to for any more complicated evaluations.

Reimplemented from NuCut.

Definition at line 2082 of file NuCutImps.cxx.

References fCC2014NewCoilCut, and NuCutImps::CC2014NewCoilCut::InFidVol().

02083   {
02084     return fCC2014NewCoilCut.InFidVol(nu);
02085   }

void NuCutImps::CC2014NewCoilCutResolution::Preselection ( const NuEvent nu  )  [virtual]

Implements NuCut.

Definition at line 2087 of file NuCutImps.cxx.

References NuCut::Defer_Preselection(), and fCC2014NewCoilCut.

02088   {
02089     Defer_Preselection(fCC2014NewCoilCut, nu);
02090   }

void NuCutImps::CC2014NewCoilCutResolution::Resolution ( const NuEvent nu  )  [protected]

Definition at line 2120 of file NuCutImps.cxx.

References NuEvent::beamTypeDB, NuEvent::energy, fQuantile, fResCuts_phe, fResCuts_r1, fResCuts_r2, fResCuts_r3, NuCut::Keep_If(), Msg::kInfo, BeamType::kL250z200i, MAXMSG, NuEvent::resolution, and NuEvent::runPeriod.

Referenced by Selection().

02121   {
02122     Double_t low_en_edge = 0;
02123     Double_t high_en_edge = 0;
02124     Double_t low_res_edge = 0;
02125     Double_t high_res_edge = 0;
02126     Double_t quantiles[6];
02127     TTree *resCuts;
02128     resCuts=fResCuts_r3;
02129     //Different cuts per run period!
02130     switch(nu.runPeriod){
02131     case 1:
02132       resCuts=fResCuts_r1;
02133       break;
02134     case 2:
02135       resCuts=fResCuts_r2;
02136       break;
02137     case 3:
02138 
02139       resCuts=fResCuts_r3;
02140       break;
02141       //default:
02142       //assert(0 && "Badly configured runPeriod");
02143     }
02144     //Whatever run period this is nowadays...
02145     if(nu.beamTypeDB == BeamType::kL250z200i) {
02146       resCuts=fResCuts_phe;
02147       MAXMSG("NuCutImps",Msg::kInfo,10)<<"Using pHE resolution cuts"<<endl;
02148     }
02149 
02150     resCuts->SetBranchAddress("quantiles", quantiles);
02151     resCuts->SetBranchAddress("lowedge", &low_en_edge);
02152     resCuts->SetBranchAddress("highedge", &high_en_edge);
02153 
02154     for(Int_t i = 0; i < resCuts->GetEntries(); i++){
02155       resCuts->GetEntry(i);
02156 
02157       if(fQuantile == 0) low_res_edge = 0.0;
02158       else low_res_edge = quantiles[fQuantile];
02159 
02160       if(fQuantile == 4) high_res_edge = 9999999.;
02161       else high_res_edge = quantiles[fQuantile+1];
02162 
02163       if(nu.energy > low_en_edge && nu.energy <= high_en_edge)
02164       {
02165         MAXMSG("NuCutImps", Msg::kInfo, 20)
02166           << "Entry:" << i << " Energy:" << nu.energy << " sigma/E:"
02167           << (nu.resolution/nu.energy) << endl
02168           << "Keep if between " << low_res_edge << " and "
02169           << high_res_edge << endl;
02170 
02171         Keep_If((nu.resolution/nu.energy) > low_res_edge &&
02172                 (nu.resolution/nu.energy) <= high_res_edge, "resolution");
02173         return;
02174        }
02175     }
02176     // God knows how we got here. We'd better keep the event, but not
02177     // double-count it, so bung it in quantile 4.
02178     Keep_If(fQuantile == 4, "Overflow");
02179   }

void NuCutImps::CC2014NewCoilCutResolution::Selection ( const NuEvent nu  )  [virtual]

Implements NuCut.

Definition at line 2092 of file NuCutImps.cxx.

References NuCut::Defer_Selection(), NuEvent::detector, fCC2014NewCoilCut, Detector::kFar, and Resolution().

02093   {
02094     Defer_Selection(fCC2014NewCoilCut, nu);
02095 
02096     if(nu.detector == Detector::kFar/* && nu.charge<0*/) Resolution(nu);
02097   }

void NuCutImps::CC2014NewCoilCutResolution::SetCutFile ( const char *  filename  )  [protected]

Definition at line 2099 of file NuCutImps.cxx.

References fResCuts_phe, fResCuts_r1, fResCuts_r2, fResCuts_r3, Msg::kInfo, Msg::kWarning, and MSG.

Referenced by CC2014NewCoilCutResolution().

02100   {
02101     TDirectory* tmpd = gDirectory; //get the right path
02102     TFile* cutfile = new TFile(filename, "READ");
02103     gDirectory = tmpd;
02104     if(cutfile->IsZombie()){
02105       MSG("NuCutImps", Msg::kWarning) << "cutfile does not exist" << endl;
02106       return;
02107     }
02108     //I know, I know, this isn't a great way to do it...
02109 
02110     fResCuts_r1 = (TTree*)cutfile->Get("cuts_r1");
02111     fResCuts_r2 = (TTree*)cutfile->Get("cuts_r2");
02112     fResCuts_r3 = (TTree*)cutfile->Get("cuts_r3");
02113     fResCuts_phe = (TTree*)cutfile->Get("cuts_phe");
02114 
02115     MSG("NuCutImps", Msg::kInfo) << "cutfile exists: "
02116                                  << fResCuts_r1->GetEntries()
02117                                  << " R1 entries" << endl;
02118   }


Member Data Documentation

Definition at line 316 of file NuCutImps.h.

Referenced by InFidVol(), Preselection(), and Selection().

Definition at line 335 of file NuCutImps.h.

Referenced by CC2014NewCoilCutResolution(), and Resolution().

Definition at line 334 of file NuCutImps.h.

Referenced by Resolution(), and SetCutFile().

Definition at line 331 of file NuCutImps.h.

Referenced by Resolution(), and SetCutFile().

Definition at line 332 of file NuCutImps.h.

Referenced by Resolution(), and SetCutFile().

Definition at line 333 of file NuCutImps.h.

Referenced by Resolution(), and SetCutFile().


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

Generated on 19 Jan 2018 for loon by  doxygen 1.6.1