NuCutImps::Preselection2014 Class Reference

#include <NuCutImps.h>

Inheritance diagram for NuCutImps::Preselection2014:
NuCut

List of all members.

Public Member Functions

 Preselection2014 (const NuPlots *plots=0)
void Preselection (const NuEvent &nu)
void Selection (const NuEvent &nu)

Private Member Functions

Bool_t InFidVol (const NuEvent &nu) const

Private Attributes

Int_t fCoilDir
 +ve -> expect !coilIsReverse. 0 -> not initialized
Int_t fBeamType

Detailed Description

Definition at line 39 of file NuCutImps.h.


Constructor & Destructor Documentation

NuCutImps::Preselection2014::Preselection2014 ( const NuPlots plots = 0  ) 

Definition at line 889 of file NuCutImps.cxx.

References NuCuts::kCC2014, NuCut::SetAnaVersion(), and NuCut::SetFidVol().

00889                                                          :
00890     NuCut("Preselection2014", plots), fCoilDir(0), fBeamType(0)
00891   {
00892     SetFidVol("cc2008");
00893     //SetAnaVersion(NuCuts::kPreselectionNC2012);
00894     SetAnaVersion(NuCuts::kCC2014);//kCC0720Std);
00895   }


Member Function Documentation

Bool_t NuCutImps::Preselection2014::InFidVol ( const NuEvent nu  )  const [private, 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 897 of file NuCutImps.cxx.

References NuEvent::anaVersion, NuEvent::detector, NuReco::GetEvtEnergy(), NuCut::InFidVolEvt(), NuLibrary::Instance(), Detector::kFar, Msg::kInfo, Detector::kNear, NuCuts::kPreselectionNC2012, Munits::m, MSG, NuEvent::nplaneShw, NuEvent::ntrk, NuLibrary::reco, NuEvent::trknplane, NuEvent::vtxMetersToCloseEdgeEvt, NuEvent::vtxMetersToCloseEdgeTrk, NuEvent::vtxMetersToCoilEvt, NuEvent::vtxMetersToCoilTrk, NuEvent::zEvtVtx, and NuEvent::zTrkVtx.

Referenced by Preselection().

00898   {
00899     //seperate bools for CCs in the CC analysis, NCs in the CC analysis and both in the NC analysis
00900     Bool_t inCC=false;
00901     Bool_t inNC=false;
00902     Bool_t inSterile=false;
00903         
00905     //cc fidvol cut version
00907     if (nu.ntrk > 0){
00908       // Move trk vtx upstream by 3.92cm from scintillator to steel
00909       // for this analysis version
00910       NuEvent nuc = nu;
00911       nuc.zTrkVtx = nu.zTrkVtx - (0.0392*Munits::m);
00912 
00913       inCC=NuCut::InFidVol(nuc);
00914     }
00915 
00917     //nc fidvol cut version
00919     // If there's a track and it's longer than the shower then use its
00920     // vertex, otherwise use the event vertex. This is what is done in
00921     // the NC analysis.
00922     if(nu.ntrk > 0 && nu.trknplane > nu.nplaneShw){
00923       // Move trk vtx upstream by 3.92cm from scintillator to steel
00924       NuEvent nuc = nu;
00925       nuc.zTrkVtx = nu.zTrkVtx - (0.0392*Munits::m);
00926       inNC=NuCut::InFidVol(nuc);
00927     }
00928     else inNC=NuCut::InFidVolEvt(nu);
00929 
00930     NuEvent nusterile = nu;
00931     nusterile.anaVersion =NuCuts::kPreselectionNC2012;
00932     NuLibrary& lib = NuLibrary::Instance();
00933     lib.reco.GetEvtEnergy(nusterile, true);
00934     
00935     //2012 NC Analysis fiducial volume cuts
00936 
00937     double distToEdge = 0.0;
00938     double distToCoil = 0.0;
00939     double vtxZ = nusterile.zEvtVtx;
00940     const int trackExtension = nusterile.nplaneShw - nusterile.trknplane;
00941 
00942     if (nusterile.ntrk > 0 && trackExtension < 0) {
00943       distToEdge = nusterile.vtxMetersToCloseEdgeTrk;
00944       distToCoil = nusterile.vtxMetersToCoilTrk;
00945       vtxZ = nusterile.zTrkVtx - 0.0392 * Munits::m;
00946     } else {
00947       distToEdge = nusterile.vtxMetersToCloseEdgeEvt;
00948       distToCoil = nusterile.vtxMetersToCoilEvt;
00949     }
00950 
00951     if (nusterile.detector == Detector::kNear) {
00952       if (distToEdge >= 0.5 && vtxZ <= 4.7368 && vtxZ >= 1.7) {
00953         inSterile = true;
00954       } 
00955     } else if (nusterile.detector == Detector::kFar) {
00956 
00957       if (distToEdge >= 0.4 && distToCoil >= 0.6) {
00958         if (vtxZ >= 0.21 && vtxZ <= 28.96) { // front / back of detector
00959           if (vtxZ <= 13.72 || vtxZ >= 16.12) { // SM gap
00960             inSterile = true;
00961           }
00962         }
00963       }
00964 
00965     } else {
00966       MSG("NuSterileCutImps", Msg::kInfo) << "Detector is not ND or FD, it is " << nusterile.detector 
00967                                           << ". No fiducial cuts applied." << endl;
00968     }
00969 
00970 
00971     if(!(inCC || inNC || inSterile)&&(inSterile)){
00972       cout<<inCC<<inNC<<inSterile<<endl;
00973       assert(0);
00974     }
00975 
00976     //return inSterile;
00977     return inCC || inNC || inSterile;
00978 
00979   }

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

Implements NuCut.

Definition at line 982 of file NuCutImps.cxx.

References bfld::AsString(), NuEvent::beamTypeDB, NuEvent::coilIsOk, NuEvent::coilIsReverse, NuCut::Cut_If(), NuEvent::cutOnDataQuality, NuEvent::cutOnSpillTiming, NuEvent::detector, fBeamType, fCoilDir, NuEvent::goodBeam, NuEvent::goodBeamSntp, NuCuts::GoodTimeToNearestSpill(), InFidVol(), NuEvent::isGoodDataQuality, NuEvent::isLI, SimFlag::kData, NuCut::Keep_Data_If(), NuCut::Keep_If(), Msg::kError, Detector::kFar, BeamType::kL010z170i, BeamType::kL010z185i, BeamType::kL010z200i, Detector::kNear, Msg::kWarning, NuEvent::litime, MAXMSG, MSG, NuEvent::run, NuEvent::simFlag, NuEvent::snarl, and NuEvent::subRun.

00983   {
00984     //The pot counting cuts on 3 variables to ensure only good spills
00985     //are included:
00986     // - beam quality
00987     // - coil quality
00988     // - detector data quality
00989     //the exact same variables have to be cut on here to be consistent
00990 
00991     Bool_t goodBeamToUse = nu.goodBeamSntp;
00992     if (Detector::kFar == nu.detector &&
00993         SimFlag::kData == nu.simFlag){
00994       goodBeamToUse = nu.goodBeam;
00995     }
00996 
00997     //only accept events where the beam hit the target well
00998     Keep_Data_If(goodBeamToUse, nu,"GoodBeam");
00999 
01000     //Is the coil okay? Only cut if we have the cutOnDataQuality flag
01001     Keep_Data_If(nu.coilIsOk || !nu.cutOnDataQuality, nu, "CoilIsOkay");
01002 
01003     //Data Quality
01004     Keep_Data_If(nu.isGoodDataQuality ||
01005                 !nu.cutOnDataQuality, nu, "DataQuality");
01006 
01007     //Filter out LI events
01008     Cut_If(nu.isLI || nu.litime != -1, "IsLI");
01009 
01010     // Per-detector cuts
01011     if (nu.detector == Detector::kFar) {
01012       //Keep events in time with the beam spill
01013       //Use the standard -2, 12us cuts
01014       Keep_Data_If(NuCuts::GoodTimeToNearestSpill(nu, -2, +12) ||
01015                   !nu.cutOnSpillTiming, nu, "SpillTime");
01016     }
01017     else if (nu.detector == Detector::kNear) {
01018       //only run the fidvol cut for the ND
01019       //since we want the RAF events at the FD
01020       Keep_If(this->InFidVol(nu),"InFidVolCCOrNC");
01021     } // End of per-detector cuts
01022 
01023     //Now do some sanity checks
01024     //files must never mix data where
01025     //- detector coil is not always the same direction
01026     //- beam type is not always the same
01027     //NOTE: for the FD data, the pot counting is less of an issue
01028     //the crucial thing is that the value of pot used was obtained
01029     //using the same cuts as are applied here
01030 
01031     if (nu.coilIsOk){
01032       //store the coil direction on the first ("good") event
01033       if (fCoilDir == 0) {
01034         if (nu.coilIsReverse) fCoilDir = -1;
01035         else fCoilDir = +1;
01036       }
01037 
01038       //check that the coil current direction is the same as first event
01039       //can't just cut it out because then the POT count will be wrong
01040       if (!((fCoilDir < 0 && nu.coilIsReverse) ||
01041             (fCoilDir > 0 && !nu.coilIsReverse)))
01042       {
01043         if (nu.detector==Detector::kNear){
01044           MAXMSG("NuCutImps",Msg::kError,100)
01045             <<"Seeing mixed coil current: fCoilDir="<<fCoilDir
01046             <<", coilIsReverse="<<nu.coilIsReverse
01047             <<", coilIsOk="<<nu.coilIsOk
01048             <<", goodBeamSntp="<<nu.goodBeamSntp
01049             <<", goodBeam="<<nu.goodBeam
01050             <<", goodBeamToUse="<<goodBeamToUse
01051             <<", isGoodDataQuality="<<nu.isGoodDataQuality
01052             <<", run="<<nu.run<<", subRun="<<nu.subRun
01053             <<", snarl="<<nu.snarl
01054             <<endl;
01055           assert((fCoilDir < 0 && nu.coilIsReverse) ||
01056                  (fCoilDir > 0 && !nu.coilIsReverse));
01057         }
01058         else if (nu.detector==Detector::kFar) {
01059           //want to print out every single one of these warnings!
01060           MSG("NuCutImps",Msg::kWarning)
01061             <<"Seeing mixed coil current: fCoilDir="<<fCoilDir
01062             <<", coilIsReverse="<<nu.coilIsReverse
01063             <<", coilIsOk="<<nu.coilIsOk
01064             <<", goodBeamSntp="<<nu.goodBeamSntp
01065             <<", goodBeam="<<nu.goodBeam
01066             <<", goodBeamToUse="<<goodBeamToUse
01067             <<", isGoodDataQuality="<<nu.isGoodDataQuality
01068             <<", run="<<nu.run<<", subRun="<<nu.subRun
01069             <<", snarl="<<nu.snarl
01070             <<endl;
01071         }
01072         else {
01073           cout<<"Ahhh, detector not known..."<<endl;
01074         }
01075       }
01076     }
01077 
01078     if (goodBeamToUse){
01079       //store the beam type on the first ("good") event
01080       if (fBeamType==0) {
01081       //treat the special horn current runs as std le-10
01082       if (nu.beamTypeDB == BeamType::kL010z170i ||
01083           nu.beamTypeDB == BeamType::kL010z200i) {
01084         fBeamType = BeamType::kL010z185i;
01085       }
01086       else fBeamType=nu.beamTypeDB;
01087     }
01088 
01089       //count the number of mixed beam types in a row
01090       static Int_t errorCounterBeamType=0;
01091       if (nu.beamTypeDB==fBeamType) errorCounterBeamType=0;
01092 
01093       //assert on the beam type being the same as the first event
01094       if (nu.beamTypeDB != fBeamType) {
01095         if ((nu.beamTypeDB == BeamType::kL010z170i ||
01096              nu.beamTypeDB == BeamType::kL010z200i)
01097             && fBeamType==BeamType::kL010z185i)
01098         {
01099             //don't assert on the special horn current runs
01100             //these beam types are ok
01101 
01102             //reset the mixed beam type counter
01103             errorCounterBeamType=0;
01104         } else {
01105           if (nu.detector==Detector::kNear){
01106             //count the number of sequential events with mixed beam type
01107             errorCounterBeamType++;
01108 
01109             //check if the number of sequential mixed beam type
01110             //events is greater than the upper limit
01111             Int_t upperLimit=10;
01112             if (errorCounterBeamType>upperLimit) {
01113               MAXMSG("NuCutImps",Msg::kError,100)
01114                 <<"Seeing mixed beam: fBeamType="<<fBeamType
01115                 << "(" << BeamType::AsString(BeamType::BeamType_t(fBeamType))
01116                 <<"), beamTypeDB="<<nu.beamTypeDB
01117                 << "(" << BeamType::AsString(BeamType::BeamType_t(nu.beamTypeDB))
01118                 <<"), goodBeamSntp="<<nu.goodBeamSntp
01119                 <<", goodBeam="<<nu.goodBeam
01120                 <<", goodBeamToUse="<<goodBeamToUse
01121                 <<", coilIsOk="<<nu.coilIsOk
01122                 <<", coilIsReverse="<<nu.coilIsReverse
01123                 <<", isGoodDataQuality="<<nu.isGoodDataQuality
01124                 <<", run="<<nu.run<<", subRun="<<nu.subRun
01125                 <<", snarl="<<nu.snarl
01126                 <<endl;
01127               MSG("NuCutImps",Msg::kError)
01128                 <<"Asserting on mixed beam "
01129                 <<" when we get more than "<<upperLimit<<" events with"
01130                 <<" mixed beam in a row"
01131                 <<", errorCounterBeamType="<<errorCounterBeamType
01132                 <<endl;
01133               assert(nu.beamTypeDB==fBeamType);
01134             }
01135             else {
01136               MSG("NuCutImps",Msg::kWarning)
01137                 <<"Not asserting on mixed beam until"
01138                 <<" we get more than "<<upperLimit<<" events with"
01139                 <<" mixed beam in a row"
01140                 <<", errorCounterBeamType="<<errorCounterBeamType
01141                 <<", fBeamType="<<fBeamType
01142                 <<"(" << BeamType::AsString(BeamType::BeamType_t(fBeamType))
01143                 <<"), beamTypeDB="<<nu.beamTypeDB
01144                 <<"(" << BeamType::AsString(BeamType::BeamType_t(nu.beamTypeDB))
01145                 <<"), goodBeamSntp="<<nu.goodBeamSntp
01146                 <<", goodBeam="<<nu.goodBeam
01147                 <<", goodBeamToUse="<<goodBeamToUse
01148                 <<", coilIsOk="<<nu.coilIsOk
01149                 <<", coilIsReverse="<<nu.coilIsReverse
01150                 <<", isGoodDataQuality="<<nu.isGoodDataQuality
01151                 <<", run="<<nu.run<<", subRun="<<nu.subRun
01152                 <<", snarl="<<nu.snarl
01153                 <<endl;
01154             }
01155           }
01156           else if (nu.detector==Detector::kFar) {
01157             //want to print out every single one of these warnings!
01158             MSG("NuCutImps",Msg::kWarning)
01159               <<"Seeing mixed beam: fBeamType="<<fBeamType
01160               <<"(" << BeamType::AsString(BeamType::BeamType_t(fBeamType))
01161               <<"), beamTypeDB="<<nu.beamTypeDB
01162               <<"(" << BeamType::AsString(BeamType::BeamType_t(nu.beamTypeDB))
01163               <<"), goodBeamSntp="<<nu.goodBeamSntp
01164               <<", goodBeam="<<nu.goodBeam
01165               <<", goodBeamToUse="<<goodBeamToUse
01166               <<", coilIsOk="<<nu.coilIsOk
01167               <<", coilIsReverse="<<nu.coilIsReverse
01168               <<", isGoodDataQuality="<<nu.isGoodDataQuality
01169               <<", run="<<nu.run<<", subRun="<<nu.subRun
01170               <<", snarl="<<nu.snarl
01171               <<endl;
01172           }
01173           else {
01174             cout<<"Ahhh, detector not known..."<<endl;
01175           }
01176         }
01177       }
01178     }
01179   }

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

Implements NuCut.

Definition at line 1182 of file NuCutImps.cxx.

01183   {
01184     // Since this is a preselection only class, we have no selection
01185   }


Member Data Documentation

Definition at line 50 of file NuCutImps.h.

Referenced by Preselection().

+ve -> expect !coilIsReverse. 0 -> not initialized

Definition at line 49 of file NuCutImps.h.

Referenced by Preselection().


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

Generated on 3 Dec 2018 for loon by  doxygen 1.6.1