MNtpModule Class Reference

#include <MNtpModule.h>

Inheritance diagram for MNtpModule:
JobCModule

List of all members.

Public Member Functions

 MNtpModule ()
 ~MNtpModule ()
void BeginJob ()
void EndJob ()
void BeginFile ()
JobCResult Ana (const MomNavigator *mom)
const RegistryDefaultConfig () const
void Config (const Registry &r)

Private Member Functions

bool PlanTriggerSim (NtpSREvent *event, NtpStRecord *record)
bool LoadLargestTrackFromEvent (NtpSREvent *event, NtpStRecord *record, int &trkidx)
bool LoadShowerAtTrackVertex (NtpSREvent *event, NtpStRecord *record, int trkidx, int &shwidx)
bool LoadLargestShowerFromEvent (NtpSREvent *event, NtpStRecord *record, int &shwidx)
void FillEvtInfo (MNtp *fMNtp, NtpSREvent *event)
void FillTrkInfo (MNtp *fMNtp, NtpSRTrack *track, NtpStRecord *record)
void FillShwInfo (MNtp *fMNtp, NtpSRShower *shower)
float CalCosine (float theta1, float phi1, float theta2, float phi2)
int Initial_state (NtpMCTruth *mctruth, TClonesArray &stdhepArray)
int Nucleus (NtpMCTruth *mctruth)
int HadronicFinalState (NtpMCTruth *mctruth, TClonesArray &stdhepArray)

Private Attributes

std::string fFile
Int_t fN
Int_t fM
bool IsTriggerOn
TFile * fNtpFile
TTree * fNtuple
Detector::Detector_t fDetector
BeamType::BeamType_t fBeam
MNtpfMNtp
MadDpID dpid
TF1 * tfit_dt_ds_pos
TF1 * tfit_dt_ds_neg
int numfiles
Detector::Detector_t fDetectorType
SimFlag::SimFlag_t fSimFlag

Detailed Description

Definition at line 34 of file MNtpModule.h.


Constructor & Destructor Documentation

MNtpModule::MNtpModule (  ) 

Definition at line 55 of file MNtpModule.cxx.

References fMNtp, tfit_dt_ds_neg, and tfit_dt_ds_pos.

00056   : fFile("mntp.root")
00057   , fN(0)
00058   , fM(0)
00059   , IsTriggerOn(false)
00060   , numfiles(0)
00061 {
00062   fMNtp = new MNtp();
00063   tfit_dt_ds_pos = new TF1("tfit_dt_ds_pos","(1./299792458.0)*x+[0]",0,30) ;
00064   tfit_dt_ds_neg = new TF1("tfit_dt_ds_neg","(-1./299792458.0)*x+[0]",0,30) ;
00065 }

MNtpModule::~MNtpModule (  ) 

Definition at line 68 of file MNtpModule.cxx.

00069 {
00070 }


Member Function Documentation

JobCResult MNtpModule::Ana ( const MomNavigator mom  )  [virtual]

Implement this for read only access to the MomNavigator

Reimplemented from JobCModule.

Definition at line 102 of file MNtpModule.cxx.

References bfld::AsString(), CalCosine(), MadDpID::CalcPID(), MadDpID::ChoosePDFs(), NtpMCGenInfo::codename, NtpSREventSummary::date, MNtp::dav_cc_pid, NtpSRDate::day, MNtp::day, NtpSRVertex::dcosx, NtpSRVertex::dcosy, NtpSRVertex::dcosz, dpid, MNtp::dt, MNtp::dt_evt, MNtp::evtazimuth, MNtp::evtcosa, NtpStRecord::evthdr, MNtp::evtno, MNtp::evtnsec, MNtp::evtsec, MNtp::evtvtxu, MNtp::evtvtxv, MNtp::evtvtxx, MNtp::evtvtxy, MNtp::evtvtxz, MNtp::evtzenith, fBeam, fDetector, fDetectorType, MNtp::fid, FillEvtInfo(), FillShwInfo(), FillTrkInfo(), NtpVtxFinder::FindVertex(), NtpMCTruth::flux, fM, fMNtp, fN, fNtuple, fSimFlag, EnergyCorrections::FullyCorrectMomentumFromRange(), EnergyCorrections::FullyCorrectShowerEnergy(), EnergyCorrections::FullyCorrectSignedMomentumFromCurvature(), NtpMCSummary::geninfo, JobCModule::Get(), VldContext::GetDetector(), SntpHelpers::GetEvent(), SntpHelpers::GetEvent2MCIndex(), MomNavigator::GetFragment(), RecRecordImp< T >::GetHeader(), MBSpill::GetHorn_current(), SntpHelpers::GetMCTruth(), MBSpill::GetNanoSec(), VldTimeStamp::GetNanoSec(), MBSpill::GetPOT(), RecDataHeader::GetRun(), RecDataHeader::GetRunType(), VldTimeStamp::GetSec(), MBSpill::GetSec(), SntpHelpers::GetShower(), VldContext::GetSimFlag(), RecPhysicsHeader::GetSnarl(), VHS::GetStrip(), RecDataHeader::GetSubRun(), MBSpill::GetTimeStamp(), VldContext::GetTimeStamp(), SntpHelpers::GetTrack(), RecPhysicsHeader::GetTrigSrc(), RecHeader::GetVldContext(), RecRecordImp< T >::GetVldContext(), NtpSRStripPulseHeight::gev, MNtp::had_fs, HadronicFinalState(), MNtp::iact, NtpMCTruth::iaction, NtpMCStdHep::IdHEP, MNtp::initial_state, Initial_state(), NtpMCTruth::inu, MNtp::inu, MNtp::ires, NtpMCTruth::iresonance, ReleaseType::IsBirch(), ReleaseType::IsCedar(), IsTriggerOn, EnergyCorrections::kBirch, CandShowerHandle::kCC, EnergyCorrections::kCedar, SimFlag::kData, Msg::kDebug, EnergyCorrections::kDefault, Msg::kError, JobCResult::kFailed, BeamType::kLE, SimFlag::kMC, Detector::kNear, JobCResult::kPassed, PlaneView::kU, PlaneView::kV, Msg::kWarning, CandShowerHandle::kWtCC, NtpSRShowerPulseHeight::linCCgev, LoadLargestShowerFromEvent(), LoadLargestTrackFromEvent(), LoadShowerAtTrackVertex(), MBSpillAccessor::LoadSpill(), AstUtil::LocalToHorizon(), ReleaseType::MakeReleaseType(), MNtp::maxuf, MNtp::maxup, MNtp::maxv, MNtp::mbhorncurr, MNtp::mbnsec, MNtp::mbpot, MNtp::mbsec, NtpStRecord::mchdr, MNtp::MN45Trig, MNtp::MNTrig, NtpSRTrack::momentum, NtpSRDate::month, MNtp::month, MSG, MNtp::muonP, NtpSREventSummary::nevent, MNtp::nevt, MNtp::nsec, Nucleus(), MNtp::nucleus, MNtp::nuenergy, MNtp::nuPx, MNtp::nuPy, MNtp::nuPz, NtpMCTruth::p4mu1, NtpMCTruth::p4neu, NtpMCTruth::p4shw, NtpMCTruth::p4tgt, MNtp::pdpx, NtpMCFluxInfo::pdpx, NtpMCFluxInfo::pdpy, MNtp::pdpy, NtpMCFluxInfo::pdpz, MNtp::pdpz, NtpSRShower::ph, phi_mb, NtpSRStrip::plane, MNtp::PlaneTrig, NtpSRStrip::planeview, PlanTriggerSim(), NtpMCFluxInfo::ptype, MNtp::ptype, MNtp::Q2, NtpMCTruth::q2, NtpSRMomentum::range, MNtp::recoEmu, MNtp::recoEshw, MNtp::recoPmu, MNtp::Reset(), MNtp::ResetAll(), MNtp::run, MNtp::runtype, MNtp::sec, EnergyCorrections::SetCorrectionVersion(), NtpVtxFinder::SetTargetEvent(), MNtp::showerE, MNtp::shwEMgev, MNtp::shwgev, MNtp::shwlinCCcor, MNtp::shwlinCCgev, MNtp::shwlinNCgev, NtpSRShower::shwph, MNtp::shwtotgev, MNtp::shwwtCCcor, MNtp::shwwtCCgev, MNtp::shwwtNCgev, MNtp::snarl, NtpStRecord::stdhep, NtpSRStrip::strip, MNtp::subrun, MNtp::tarE, MNtp::tarPx, MNtp::tarPy, MNtp::tarPz, theta_mb, NtpStRecord::thtrk, NtpSRStrip::time1, MNtp::tptype, NtpMCFluxInfo::tptype, MNtp::tpx, NtpMCFluxInfo::tpx, NtpMCFluxInfo::tpy, MNtp::tpy, NtpMCFluxInfo::tpz, MNtp::tpz, MNtp::trgsrc, MNtp::trkazimuth, MNtp::trkeqp, MNtp::trkexit, MNtp::trkp_fit, MNtp::trkp_fitcor, MNtp::trkp_range, MNtp::trkp_rangecor, MNtp::trkpassfit, MNtp::trkpid, MNtp::trkqp, NtpTHTrack::trkstdhep, MNtp::trkzenith, MNtp::tvx, NtpMCFluxInfo::tvx, MNtp::tvy, NtpMCFluxInfo::tvy, MNtp::tvz, NtpMCFluxInfo::tvz, NtpSRShower::vtx, NtpSRTrack::vtx, NtpVtxFinder::VtxU(), NtpVtxFinder::VtxV(), NtpVtxFinder::VtxX(), NtpVtxFinder::VtxY(), NtpVtxFinder::VtxZ(), MNtp::vx, NtpMCFluxInfo::vx, MNtp::vy, NtpMCFluxInfo::vy, MNtp::vz, NtpMCFluxInfo::vz, NtpMCTruth::w2, MNtp::W2, NtpMCTruth::x, MNtp::x, NtpMCTruth::y, and MNtp::y.

00103 {
00104   //reset everything at the beginning of each snarl
00105   fMNtp->ResetAll();
00106 
00107   //get run no,subrun no and snarl no
00108   RecRecordImp<RecCandHeader> *rr = 
00109     dynamic_cast<RecRecordImp<RecCandHeader>*>
00110     (mom->GetFragment("RecRecordImp<RecCandHeader>"));
00111   if (rr) {
00112     fMNtp->run    = rr->GetHeader().GetRun();
00113     fMNtp->subrun = rr->GetHeader().GetSubRun();
00114     fMNtp->snarl  = rr->GetHeader().GetSnarl();
00115     fDetectorType = rr->GetHeader().GetVldContext().GetDetector();
00116     fSimFlag = rr->GetHeader().GetVldContext().GetSimFlag();
00117   }
00118   else {
00119     MSG("MNtpModule", Msg::kWarning) << "Could not get RecCandHeader from MOM" << endl;
00120   }
00121 
00122 
00123   //get NtpStRecord
00124   NtpStRecord *record = dynamic_cast<NtpStRecord *>(mom->GetFragment("NtpStRecord"));
00125   if (!record) {
00126     MSG("MNtpModule", Msg::kError) << "Could not get NtpStRecord from MOM" << endl;
00127     return JobCResult::kFailed;
00128   }
00129 
00130   //get no. of events
00131   int nevt = record->evthdr.nevent;
00132   if (!nevt) {
00133     MSG("MNtpModule", Msg::kDebug) << "There is no event in this snarl"<<endl;
00134     return JobCResult::kPassed;
00135   }
00136 
00137   fMNtp->nevt = nevt;
00138 
00139   //get VldContext, UgliGeomHandle and VldTimeStamp
00140   VldContext vldc = *(record->GetVldContext());
00141   VldTimeStamp vldts = vldc.GetTimeStamp();
00142   //UgliGeomHandle ugh(vldc);
00143 
00144   fDetector = record->GetHeader().GetVldContext().GetDetector();
00145   fBeam = BeamType::kLE ;
00146 
00147   fMNtp->sec = (int)vldts.GetSec();
00148   fMNtp->nsec = (int)vldts.GetNanoSec();
00149 
00150   fMNtp->day = (int)record->evthdr.date.day;
00151   fMNtp->month = (int)record->evthdr.date.month;
00152 
00153   //get trigger source:
00154   fMNtp->trgsrc = record->GetHeader().GetTrigSrc();
00155 
00156   //get run type:
00157   fMNtp->runtype = record->GetHeader().GetRunType();
00158 
00159   string relName = record->GetTitle();
00160   string mcinfo = "";
00161   if (fSimFlag==SimFlag::kMC){
00162     mcinfo = "Carrot";
00163     string temp = record->mchdr.geninfo.codename;
00164     if (temp.size() != 0) mcinfo = temp;
00165     if (mcinfo == "development") mcinfo = "daikon_04";
00166   }
00167   int release = ReleaseType::MakeReleaseType(relName, mcinfo);
00168   string title = ReleaseType::AsString(release);
00169   static bool checkrel = true;
00170   if (checkrel){
00171     cout<<"Release "<<release<<": "<<title<<endl;
00172     checkrel = false;
00173   }
00174   // ENERGY CORRECTIONS - 1) set reco version
00175   if(ReleaseType::IsCedar(release))
00176     EnergyCorrections::SetCorrectionVersion(EnergyCorrections::kCedar);
00177   if(ReleaseType::IsBirch(release))
00178     EnergyCorrections::SetCorrectionVersion(EnergyCorrections::kBirch);
00179 
00180   // ENERGY CORRECTIONS - 2) set release type for MC or data
00181 
00182   // MC or data
00183   ReleaseType::Release_t reltype=release;
00184 
00185 
00186   // ENERGY CORRECTIONS - 3) set correction version - use default
00187 
00188   EnergyCorrections::WhichCorrection_t corrver=EnergyCorrections::kDefault;
00189 
00190   
00191   
00192   if (fSimFlag==SimFlag::kData && !IsTriggerOn){
00193     //get MinoBooNE spill
00194     MBSpillAccessor& spillAccess = MBSpillAccessor::Get();
00195     
00196     const MBSpill *spillInfo = spillAccess.LoadSpill(vldc.GetTimeStamp());
00197     if (spillInfo) {
00198       fMNtp->dt = double(spillInfo->GetTimeStamp()-vldc.GetTimeStamp());
00199       fMNtp->mbsec = spillInfo->GetSec();
00200       fMNtp->mbnsec = spillInfo->GetNanoSec();
00201       fMNtp->mbpot = spillInfo->GetPOT();
00202       fMNtp->mbhorncurr = spillInfo->GetHorn_current();
00203     }
00204   }
00205 
00206   for(int ievt = 0; ievt<nevt; ievt++){//loop over all events
00207     fMNtp->Reset();  //reset all the variables excepts the header information
00208     fMNtp->evtno = ievt;
00209     NtpSREvent *event = SntpHelpers::GetEvent(ievt,record);
00210     if (!event) continue;
00211 
00212     if (IsTriggerOn){
00213       fN = 10 ;
00214       fM = 12 ;
00215       fMNtp->MNTrig = PlanTriggerSim(event,record) ;
00216 
00217       fN = 10 ;
00218       fM = 600 ;
00219       fMNtp->PlaneTrig = PlanTriggerSim(event,record) ;
00220 
00221       fN = 4 ;
00222       fM = 5 ;
00223       fMNtp->MN45Trig = PlanTriggerSim(event,record) ; 
00224     }
00225 
00226     //index will be -1 if no track/shower present in event
00227     //if there is a track, find the primary shower close to the track vertex
00228     //if there is no track, find the biggest shower
00229     //pull out from Mad
00230     int track_index   = -1;
00231     int shower_index  = -1;
00232     if(ReleaseType::IsBirch(release)){
00233       if (LoadLargestTrackFromEvent(event,record,track_index)){
00234         if(!LoadShowerAtTrackVertex(event,record,track_index,shower_index)){
00235           LoadLargestShowerFromEvent(event,record,shower_index);
00236         }
00237       }
00238       else LoadLargestShowerFromEvent(event,record,shower_index);
00239     }
00240     else if(ReleaseType::IsCedar(release)){
00241       if (event->ntrack){
00242         track_index = event->trk[0];
00243       }
00244       if (event->nshower){
00245         shower_index = event->shw[0];
00246       }
00247     }
00248 
00249     if (track_index!=-1){//found track
00250       NtpSRTrack *track = SntpHelpers::GetTrack(track_index,record);
00251       if (track) {
00252         this->FillTrkInfo(fMNtp,track,record);
00253         fMNtp->trkp_rangecor = FullyCorrectMomentumFromRange(fMNtp->trkp_range,vldc,reltype,corrver);
00254         fMNtp->trkp_fitcor = FullyCorrectSignedMomentumFromCurvature(fMNtp->trkp_fit,vldc,reltype,corrver);      
00255         fMNtp->recoPmu = 0.;
00256         if (fMNtp->trkp_rangecor>0){
00257           fMNtp->recoPmu = fMNtp->trkp_rangecor;
00258         }
00259         if (fMNtp->trkpassfit&&fMNtp->trkexit&&TMath::Abs(fMNtp->trkeqp/fMNtp->trkqp)<0.3){
00260           fMNtp->recoPmu = TMath::Abs(fMNtp->trkp_fitcor);
00261         }
00262         
00263         fMNtp->recoEmu = sqrt(pow(fMNtp->recoPmu,2)+0.105658*0.105658);
00264       }
00265     }
00266 
00267     if (shower_index!=-1){//found shower
00268       NtpSRShower *shower = SntpHelpers::GetShower(shower_index,record);
00269       if (shower) {
00270         this->FillShwInfo(fMNtp,shower);
00271         fMNtp->shwlinCCcor = FullyCorrectShowerEnergy(fMNtp->shwlinCCgev,CandShowerHandle::kCC,vldc,reltype,corrver);
00272         fMNtp->shwwtCCcor = FullyCorrectShowerEnergy(fMNtp->shwwtCCgev,CandShowerHandle::kWtCC,vldc,reltype,corrver);
00273         fMNtp->recoEshw = 0;
00274         if (fMNtp->shwlinCCcor>0) fMNtp->recoEshw = fMNtp->shwlinCCcor;
00275       }
00276     }
00277     else {
00278       fMNtp->shwgev = 0;
00279       fMNtp->shwlinCCgev = 0;
00280       fMNtp->shwwtCCgev = 0;
00281       fMNtp->shwlinCCcor = 0;
00282       fMNtp->shwwtCCcor = 0;
00283       fMNtp->shwlinNCgev = 0;
00284       fMNtp->shwwtNCgev = 0;
00285       fMNtp->shwEMgev = 0;
00286       fMNtp->shwtotgev = 0;
00287       fMNtp->recoEshw = 0;
00288     }
00289 
00290     //calculate david's cc pid
00291     if (track_index!=-1){
00292       NtpSRTrack *track = SntpHelpers::GetTrack(track_index,record);
00293       if(dpid.ChoosePDFs(fDetector,fBeam))  // david pid 
00294          fMNtp->dav_cc_pid = dpid.CalcPID(track,event,
00295                                   &(record->evthdr),fDetector,0);
00296     }  
00297 
00298     //fill event information
00299     this->FillEvtInfo(fMNtp,event);
00300     //correcting vtx
00301     if(ReleaseType::IsCedar(release)){
00302       NtpVtxFinder vtxf;
00303       vtxf.SetTargetEvent(ievt,record);
00304       if (vtxf.FindVertex()>0){
00305         fMNtp->evtvtxx = vtxf.VtxX();
00306         fMNtp->evtvtxy = vtxf.VtxY();
00307         fMNtp->evtvtxz = vtxf.VtxZ();
00308         fMNtp->evtvtxu = vtxf.VtxU();
00309         fMNtp->evtvtxv = vtxf.VtxV();
00310         if (fMNtp->evtvtxz>0.5&&fMNtp->evtvtxz<6&&
00311             pow(fMNtp->evtvtxx-1.4885,2)+pow(fMNtp->evtvtxy-0.1397,2)<1){
00312           fMNtp->fid = 1;
00313         }//fiducial volume cut
00314         else fMNtp->fid = 0;
00315       }
00316     }
00317 
00318 
00319    //calculate event direction: track+shower
00320     if (track_index!=-1&&shower_index!=-1){
00321       NtpSRTrack *track = SntpHelpers::GetTrack(track_index,record);
00322       NtpSRShower *shower = SntpHelpers::GetShower(shower_index,record);
00323       float trkdx = track->vtx.dcosx;
00324       float trkdy = track->vtx.dcosy;
00325       float trkdz = track->vtx.dcosz;
00326       float shwdx = shower->vtx.dcosx;
00327       float shwdy = shower->vtx.dcosy;
00328       float shwdz = shower->vtx.dcosz;
00329       float trkmom = track->momentum.range;
00330       float shwmom = shower->ph.gev;
00331       if (shower->shwph.linCCgev>0) shwmom = shower->shwph.linCCgev;
00332       float px = trkmom*trkdx + shwmom*shwdx;
00333       float py = trkmom*trkdy + shwmom*shwdy;
00334       float pz = trkmom*trkdz + shwmom*shwdz;      
00335 
00336       if (px||py||pz){
00337         float norm = sqrt(px*px+py*py+pz*pz);
00338         float dcosx = -px/norm;
00339         float dcosy = -py/norm;
00340         float dcosz = -pz/norm;
00341         double altitude, evtazimuth_tmp;
00342         AstUtil::LocalToHorizon(dcosx,dcosy,dcosz,Detector::kNear,altitude,evtazimuth_tmp);
00343         fMNtp->evtzenith = 90. - altitude;
00344         fMNtp->evtazimuth = evtazimuth_tmp;
00345         fMNtp->evtcosa = CalCosine(fMNtp->evtzenith,fMNtp->evtazimuth,theta_mb,phi_mb);
00346       }
00347     }
00348     else if (track_index!=-1){
00349       fMNtp->evtzenith = fMNtp->trkzenith;
00350       fMNtp->evtazimuth = fMNtp->trkazimuth;
00351       fMNtp->evtcosa = CalCosine(fMNtp->evtzenith,fMNtp->evtazimuth,theta_mb,phi_mb);
00352 
00353     }
00354     
00355     int nstp = event->nstrip; //no. of strips
00356       
00357     // the following is to check the timing of the earliest hit of an event is consistent with
00358     // the snarl time (just for MB events)
00359     if (fSimFlag==SimFlag::kData && !IsTriggerOn){
00360       //find the earliest time in the event
00361       //because Peter worried there might be a difference between 
00362       //the snarl timestamp and the earliest hit time.
00363       VldTimeStamp t0(2020,1,1,0,0,0);
00364       for (int istp = 0; istp<nstp; istp++){
00365         NtpSRStrip *strip = SntpHelpers::GetStrip(event->stp[istp],record);
00366         //cout<<strip->time1<<" "<<VldTimeStamp(strip->time1)<<" "<<t0<<endl;
00367         if (fMNtp->sec+strip->time1<double(t0)){
00368           t0 = VldTimeStamp(fMNtp->sec+strip->time1);
00369         }
00370       }
00371       
00372       fMNtp->evtsec = (int)t0.GetSec();
00373       fMNtp->evtnsec = (int)t0.GetNanoSec();
00374 
00375       MBSpillAccessor& spillAccess = MBSpillAccessor::Get();
00376 
00377       const MBSpill *spillInfo0 = spillAccess.LoadSpill(t0);
00378       if (spillInfo0) {
00379         fMNtp->dt_evt = double(spillInfo0->GetTimeStamp()-t0);
00380       }
00381 //      if (spillInfo&&spillInfo0&&spillInfo->GetTimeStamp()!=spillInfo0->GetTimeStamp()){
00382 //      MSG("MNtpModule", Msg::kWarning) <<"Different time stamps between snarl and earliest hit "<<vldts<<" "<<t0<<endl;
00383 //      }
00384     }
00385   
00386     //find the highest strip number in fu,pu and v views.
00387     //this is not as useful as trkendu,trkendv and trkendy
00388     for (int istp = 0; istp < nstp; istp++){
00389       NtpSRStrip *strip = SntpHelpers::GetStrip(event->stp[istp],record);
00390       if(strip->planeview==PlaneView::kU){
00391         if((strip->plane-1)%5==0){//fully instrumented
00392           if (strip->strip>fMNtp->maxuf) fMNtp->maxuf = strip->strip;
00393         }
00394         else {//partially instrumented
00395           if (strip->strip>fMNtp->maxup) fMNtp->maxup = strip->strip;
00396         }
00397       }
00398       else if (strip->planeview==PlaneView::kV){
00399         if (strip->strip>fMNtp->maxv) fMNtp->maxv = strip->strip;
00400       }
00401     }
00402     
00403     //fill mc information
00404     if (fSimFlag==SimFlag::kMC){//mc
00405       int index = SntpHelpers::GetEvent2MCIndex(fMNtp->evtno,record);
00406       NtpMCTruth *mctruth = SntpHelpers::GetMCTruth(index,record);
00407       TClonesArray& stdhepArray = *(record->stdhep);
00408       if (mctruth) {
00409         fMNtp->nuenergy = mctruth->p4neu[3];
00410         fMNtp->nuPx = mctruth->p4neu[0];
00411         fMNtp->nuPy = mctruth->p4neu[1];
00412         fMNtp->nuPz = mctruth->p4neu[2];
00413         fMNtp->tarE = mctruth->p4tgt[3];
00414         fMNtp->tarPx = mctruth->p4tgt[0];
00415         fMNtp->tarPy = mctruth->p4tgt[1];
00416         fMNtp->tarPz = mctruth->p4tgt[2];
00417         fMNtp->ires = mctruth->iresonance;
00418         fMNtp->inu = mctruth->inu;
00419         fMNtp->iact = mctruth->iaction;
00420         fMNtp->initial_state = this->Initial_state(mctruth,stdhepArray);
00421         fMNtp->nucleus = this->Nucleus(mctruth);
00422         fMNtp->had_fs = this->HadronicFinalState(mctruth,stdhepArray);
00423         fMNtp->ptype = mctruth->flux.ptype;
00424         fMNtp->tptype = mctruth->flux.tptype;
00425         fMNtp->pdpx = mctruth->flux.pdpx;
00426         fMNtp->pdpy = mctruth->flux.pdpy;
00427         fMNtp->pdpz = mctruth->flux.pdpz;
00428         fMNtp->vx = mctruth->flux.vx;
00429         fMNtp->vy = mctruth->flux.vy;
00430         fMNtp->vz = mctruth->flux.vz;
00431         fMNtp->tpx = mctruth->flux.tpx;
00432         fMNtp->tpy = mctruth->flux.tpy;
00433         fMNtp->tpz = mctruth->flux.tpz;
00434         fMNtp->tvx = mctruth->flux.tvx;
00435         fMNtp->tvy = mctruth->flux.tvy;
00436         fMNtp->tvz = mctruth->flux.tvz;
00437         fMNtp->x = mctruth->x;
00438         fMNtp->y = mctruth->y;
00439         fMNtp->Q2 = mctruth->q2;
00440         fMNtp->W2 = mctruth->w2;
00441         fMNtp->showerE = mctruth->p4shw[3];
00442         if (pow(mctruth->p4mu1[0],2)+pow(mctruth->p4mu1[1],2)+pow(mctruth->p4mu1[2],2)>0) fMNtp->muonP = sqrt(pow(mctruth->p4mu1[0],2)+pow(mctruth->p4mu1[1],2)+pow(mctruth->p4mu1[2],2));
00443       }
00444       if (track_index!=-1){
00445         NtpTHTrack *thtrack = dynamic_cast<NtpTHTrack *>((*record->thtrk)[track_index]);
00446         if (thtrack){
00447           NtpMCStdHep *stdhep = dynamic_cast<NtpMCStdHep *>((*record->stdhep)[thtrack->trkstdhep]);
00448           if (stdhep){
00449             fMNtp->trkpid = stdhep->IdHEP;
00450           }
00451         }
00452       }
00453     }    
00454 
00455     fNtuple->Fill();
00456 
00457   }
00458   return JobCResult::kPassed; // kNoDecision, kFailed, etc.
00459 }

void MNtpModule::BeginFile (  )  [virtual]

Implement for notification of begin of file. See GetCurrentFile().

Reimplemented from JobCModule.

Definition at line 93 of file MNtpModule.cxx.

References numfiles.

00094 {
00095   numfiles++;
00096   cout <<"=====================================================" <<endl;
00097   cout << "This is the NO. " << numfiles << " ntuples." << endl;  
00098 }

void MNtpModule::BeginJob ( void   )  [virtual]

Implement for notification of begin of job

Reimplemented from JobCModule.

Definition at line 74 of file MNtpModule.cxx.

References fFile, fMNtp, fNtpFile, and fNtuple.

00075 {
00076   fNtpFile = new TFile(fFile.c_str(),"RECREATE");  
00077   fNtuple = new TTree("mini","Analysis Tree");
00078   fNtuple->Branch("mntp","MNtp",&fMNtp,64000,2);
00079 
00080 
00081 }

float MNtpModule::CalCosine ( float  theta1,
float  phi1,
float  theta2,
float  phi2 
) [private]

Definition at line 1173 of file MNtpModule.cxx.

Referenced by Ana(), and FillTrkInfo().

01174 {
01175   float t1 = theta1*3.1416/180;
01176   float t2 = theta2*3.1416/180;
01177   float p1 = phi1*3.1416/180;
01178   float p2 = phi2*3.1416/180;
01179   return TMath::Sin(t1)*TMath::Sin(t2)*TMath::Cos(p1-p2)+TMath::Cos(t1)*TMath::Cos(t2);
01180 }

void MNtpModule::Config ( const Registry r  )  [virtual]

Return the actual configuration. If your module directly pulls its configuration from the fConfig Registry, you don't need to override this. Override if you have local config variables.

Configure the module given the Registry r

Reimplemented from JobCModule.

Definition at line 1205 of file MNtpModule.cxx.

References fFile, Registry::GetCharString(), Registry::GetInt(), and IsTriggerOn.

01206 {
01210   fFile = r.GetCharString("File");
01211   IsTriggerOn = (bool)r.GetInt("TriggerOn") ;
01212 }

const Registry & MNtpModule::DefaultConfig ( void   )  const [virtual]

Get the default configuration registry. This should normally be overridden. One useful idiom is to implement it like:

const Registry& MyModule::DefaultConfig() const { static Registry cfg; // never is destroyed if (cfg.Size()) return cfg; // already filled it // set defaults: cfg.Set("TheAnswer",42); cfg.Set("Units","unknown"); return cfg; }

Supply the default configuration for the module

Reimplemented from JobCModule.

Definition at line 1182 of file MNtpModule.cxx.

References JobCModule::GetName(), Registry::LockValues(), Registry::Set(), and Registry::UnLockValues().

01183 {
01187   static Registry r; // Default configuration for module
01188 
01189   // Set name of config
01190   std::string name = this->GetName();
01191   name += ".config.default";
01192   r.SetName(name.c_str());
01193 
01194   // Set values in configuration
01195   r.UnLockValues();
01196   r.Set("File","mntp.root");
01197   r.Set("TriggerOn",0) ;
01198   r.LockValues();
01199 
01200   return r;
01201 }

void MNtpModule::EndJob (  )  [virtual]

Implement for notification of end of job

Reimplemented from JobCModule.

Definition at line 85 of file MNtpModule.cxx.

References fNtpFile, and fNtuple.

00086 {
00087   fNtuple->Write();
00088   fNtpFile->Close();
00089 }

void MNtpModule::FillEvtInfo ( MNtp fMNtp,
NtpSREvent event 
) [private]

Definition at line 1033 of file MNtpModule.cxx.

References MNtp::evtvtxu, MNtp::evtvtxv, MNtp::evtvtxx, MNtp::evtvtxy, MNtp::evtvtxz, MNtp::fid, MNtp::nshw, MNtp::ntrk, MNtp::recoEmu, MNtp::recoEnu, MNtp::recoEshw, MNtp::recoPmu, MNtp::recoQ2, MNtp::recoW2, MNtp::trkcosa, NtpSREvent::vtx, NtpSRVertex::x, NtpSRVertex::y, and NtpSRVertex::z.

Referenced by Ana().

01034 {
01035   fMNtp->ntrk = event->ntrack;
01036   fMNtp->nshw = event->nshower;
01037   fMNtp->evtvtxx = event->vtx.x;
01038   fMNtp->evtvtxy = event->vtx.y;
01039   fMNtp->evtvtxz = event->vtx.z;
01040   fMNtp->evtvtxu = event->vtx.u;
01041   fMNtp->evtvtxv = event->vtx.v;
01042 
01043   if (event->vtx.z>0.5&&event->vtx.z<6&&
01044       pow(event->vtx.x-1.4885,2)+pow(event->vtx.y-0.1397,2)<1){
01045     fMNtp->fid = 1;
01046   }//fiducial volume cut
01047   else fMNtp->fid = 0;
01048 
01049   fMNtp->recoEnu = fMNtp->recoEmu;
01050   if (fMNtp->recoEshw>0) fMNtp->recoEnu += fMNtp->recoEshw;
01051 
01052   if (fMNtp->recoPmu>0){
01053     fMNtp->recoQ2 = 2*fMNtp->recoEnu*fMNtp->recoEmu-2*fMNtp->recoPmu*fMNtp->recoEnu*fMNtp->trkcosa - 0.105658*0.105658;
01054     float M = (0.9383+0.9396)/2;
01055     fMNtp->recoW2 = M*M+2*M*fMNtp->recoEshw-fMNtp->recoQ2;
01056   }
01057 }

void MNtpModule::FillShwInfo ( MNtp fMNtp,
NtpSRShower shower 
) [private]

Definition at line 1152 of file MNtpModule.cxx.

References NtpSRShowerPulseHeight::EMgev, NtpSRStripPulseHeight::gev, NtpSRShowerPulseHeight::linCCgev, NtpSRShowerPulseHeight::linNCgev, NtpSRShower::ph, MNtp::shwEMgev, MNtp::shwgev, MNtp::shwlinCCgev, MNtp::shwlinNCgev, NtpSRShower::shwph, MNtp::shwvtxu, MNtp::shwvtxv, MNtp::shwvtxx, MNtp::shwvtxy, MNtp::shwvtxz, MNtp::shwwtCCgev, MNtp::shwwtNCgev, NtpSRVertex::u, NtpSRVertex::v, NtpSRShower::vtx, NtpSRShowerPulseHeight::wtCCgev, NtpSRShowerPulseHeight::wtNCgev, NtpSRVertex::x, NtpSRVertex::y, and NtpSRVertex::z.

Referenced by Ana().

01153 {
01154   fMNtp->shwvtxx = shower->vtx.x;
01155   fMNtp->shwvtxy = shower->vtx.y;
01156   fMNtp->shwvtxz = shower->vtx.z;
01157   fMNtp->shwvtxu = shower->vtx.u;
01158   fMNtp->shwvtxv = shower->vtx.v;
01159   fMNtp->shwgev = shower->ph.gev;
01160   fMNtp->shwlinCCgev = shower->shwph.linCCgev;
01161   fMNtp->shwwtCCgev = shower->shwph.wtCCgev;
01162   fMNtp->shwlinNCgev = shower->shwph.linNCgev;
01163   fMNtp->shwwtNCgev = shower->shwph.wtNCgev;
01164   fMNtp->shwEMgev = shower->shwph.EMgev;
01165   //make corrections for shower energy
01166 //  fMNtp->shwlinCCcor = CorrectShowerEnergy(fMNtp->shwlinCCgev,fDetectorType,CandShowerHandle::kCC); 
01167 //  fMNtp->shwwtCCcor = CorrectShowerEnergy(fMNtp->shwwtCCgev,fDetectorType,CandShowerHandle::kWtCC); 
01168 
01169 }

void MNtpModule::FillTrkInfo ( MNtp fMNtp,
NtpSRTrack track,
NtpStRecord record 
) [private]

Definition at line 1061 of file MNtpModule.cxx.

References NtpSRCosmicRay::azimuth, NtpSRPlane::beg, CalCosine(), UgliStripHandle::ClearFiber(), NtpSRTrack::cr, NtpSRFiducial::dr, NtpSRTrack::ds, NtpSRPlane::end, NtpSRTrack::end, NtpSRMomentum::eqp, NtpSRTrack::fidall, NtpSRTrack::fit, UgliStripHandle::GetHalfLength(), VHS::GetStrip(), UgliGeomHandle::GetStripHandle(), RecRecordImp< T >::GetVldContext(), UgliStripHandle::GlobalToLocal(), Detector::kNear, StripEnd::kWest, NtpSRTrack::momentum, NtpSRPlane::n, NtpSRTrack::nstrip, NtpSRFitTrack::pass, phi_mb, NtpSRStrip::plane, NtpSRTrack::plane, NtpSRMomentum::qp, NtpSRMomentum::range, NtpSRTrack::stp, NtpSRTrack::stpds, NtpSRTrack::stpx, NtpSRTrack::stpy, NtpSRTrack::stpz, NtpSRStrip::strip, tfit_dt_ds_neg, tfit_dt_ds_pos, theta_mb, NtpSRStrip::time1, MNtp::trkazimuth, MNtp::trkcosa, MNtp::trkendu, MNtp::trkendv, MNtp::trkendx, MNtp::trkendy, MNtp::trkendz, MNtp::trkeqp, MNtp::trkexit, MNtp::trkfiddr, MNtp::trklength, MNtp::trkp_fit, MNtp::trkp_range, MNtp::trkpassfit, MNtp::trkplanes, MNtp::trkqp, MNtp::trktrms_neg, MNtp::trktrms_pos, MNtp::trkvtxu, MNtp::trkvtxv, MNtp::trkvtxx, MNtp::trkvtxy, MNtp::trkvtxz, MNtp::trkzenith, NtpSRVertex::u, NtpSRVertex::v, PropagationVelocity::Velocity(), NtpSRTrack::vtx, UgliStripHandle::WlsPigtail(), NtpSRVertex::x, NtpSRVertex::y, NtpSRVertex::z, and NtpSRCosmicRay::zenith.

Referenced by Ana().

01062 {
01063   fMNtp->trkplanes = track->plane.n;
01064   fMNtp->trklength = abs(track->plane.end-track->plane.beg)+1;
01065   fMNtp->trkp_range = track->momentum.range;
01066 //  bool isdata=false;
01067 //  if (fSimFlag == SimFlag::kData) isdata=true;
01068 //  fMNtp->trkp_rangecor = CorrectMomentumFromRange(fMNtp->trkp_range, isdata, fDetectorType);
01069   if (track->fit.pass){
01070     fMNtp->trkpassfit = 1;
01071   }
01072   else fMNtp->trkpassfit = 0;
01073   fMNtp->trkqp = track->momentum.qp;
01074   fMNtp->trkeqp = track->momentum.eqp;
01075   if (fMNtp->trkqp){
01076     fMNtp->trkp_fit = 1/fMNtp->trkqp;
01077     //fMNtp->trkp_fitcor = CorrectSignedMomentumFromCurvature(fMNtp->trkp_fit,isdata, fDetectorType);
01078   }
01079   fMNtp->trkvtxx = track->vtx.x;
01080   fMNtp->trkvtxy = track->vtx.y;
01081   fMNtp->trkvtxz = track->vtx.z;
01082   fMNtp->trkvtxu = track->vtx.u;
01083   fMNtp->trkvtxv = track->vtx.v;
01084   fMNtp->trkendx = track->end.x;
01085   fMNtp->trkendy = track->end.y;
01086   fMNtp->trkendz = track->end.z;
01087   fMNtp->trkendu = track->end.u;
01088   fMNtp->trkendv = track->end.v;
01089   fMNtp->trkzenith = track->cr.zenith;
01090   fMNtp->trkazimuth = track->cr.azimuth;
01091   fMNtp->trkcosa = CalCosine(fMNtp->trkzenith,fMNtp->trkazimuth,theta_mb,phi_mb);
01092   fMNtp->trkfiddr = track->fidall.dr;
01093   //determine wheather the track exits the detector
01094   //doc-1547 p6
01095   bool out = true;
01096   //"pitt" fiducial volume
01097   if (track->end.u>0.3&&track->end.u<1.8
01098       &&track->end.v>-1.8&&track->end.v<-0.3&&track->end.x<2.4
01099       &&pow(track->end.x,2)+pow(track->end.y,2)>0.64) out = false;
01100   if ((track->end.z < 7&&out)
01101       ||(track->end.z > 7&&pow(track->end.x-0.8,2)/(1.7*1.7)+pow(track->end.y,2)/(1.4*1.4)>1&&pow(track->end.x,2)+pow(track->end.y,2)>1)
01102       ||(track->end.z>15.6)){
01103     fMNtp->trkexit = 1;
01104   }
01105   else fMNtp->trkexit = 0;
01106   //cout<<out<<" "<<track->end.z<<" "<<pow(track->end.x-0.8,2)/(1.7*1.7)+pow(track->end.y,2)/(1.4*1.4)<<" "<<pow(track->end.x,2)+pow(track->end.y,2)<<" "<<fMNtp->trkexit<<endl;
01107 
01108   //calculate rms from dt_ds fitting
01109   //get VldContext, UgliGeomHandle and VldTimeStamp
01110   VldContext vldc = *(record->GetVldContext());
01111   UgliGeomHandle ugh(vldc);  vector<double> spathLength;
01112   vector<double> st0;
01113   int nstp = track->nstrip;
01114   for (int i = 0; i<nstp; i++){
01115     spathLength.push_back(track->ds-track->stpds[i]);
01116     NtpSRStrip *strip = SntpHelpers::GetStrip(track->stp[i],record);
01117     PlexStripEndId seid(Detector::kNear,strip->plane,strip->strip,StripEnd::kWest);
01118     UgliStripHandle stripHandle = ugh.GetStripHandle(seid);
01119     float halfLength = stripHandle.GetHalfLength();
01120     const TVector3 ghitxyz(track->stpx[i],track->stpy[i],track->stpz[i]);
01121     TVector3 lhitxyz = stripHandle.GlobalToLocal(ghitxyz);
01122     float fiberDist = (halfLength - lhitxyz.x() + stripHandle.ClearFiber(StripEnd::kWest) + stripHandle.WlsPigtail(StripEnd::kWest));
01123     //st0.push_back(track->stpt1[i]-fiberDist/PropagationVelocity::Velocity());
01124     st0.push_back(strip->time1-fiberDist/PropagationVelocity::Velocity());
01125   }
01126   TGraph *gr = new TGraph(st0.size(),&spathLength[0],&st0[0]);
01127   tfit_dt_ds_neg -> SetParameter(0,0.0) ;
01128   gr->Fit("tfit_dt_ds_neg","QR") ;
01129   
01130   tfit_dt_ds_pos -> SetParameter(0,0.0) ;
01131   gr->Fit("tfit_dt_ds_pos","QR") ;
01132   
01133   double trms1 = 0;
01134   double trms2 = 0;
01135   
01136   for (unsigned i = 0; i<st0.size(); i++){
01137     double x,y;
01138     gr->GetPoint(i,x,y);
01139     trms1 += pow(y-tfit_dt_ds_pos->Eval(x),2);
01140     trms2 += pow(y-tfit_dt_ds_neg->Eval(x),2);
01141   }
01142   trms1/=st0.size();
01143   trms2/=st0.size();
01144   trms1=sqrt(trms1)*1e9;
01145   trms2=sqrt(trms2)*1e9;
01146   fMNtp->trktrms_pos = trms1;
01147   fMNtp->trktrms_neg = trms2;
01148 }

int MNtpModule::HadronicFinalState ( NtpMCTruth mctruth,
TClonesArray &  stdhepArray 
) [private]

Definition at line 861 of file MNtpModule.cxx.

References NtpMCStdHep::IdHEP, NtpMCTruth::index, NtpMCTruth::iresonance, NtpMCStdHep::IstHEP, and NtpMCStdHep::mc.

Referenced by Ana().

00862                                                              {
00863   Int_t hfs=0;
00864   Int_t proc = mctruth->iresonance;
00865   Int_t nStdHep = stdhepArray.GetEntries();
00866   if(proc==1002){
00867     for(int i=0;i<nStdHep;i++){
00868 
00869       NtpMCStdHep *ntpStdHep = dynamic_cast<NtpMCStdHep *>(stdhepArray[i]);
00870 
00871       if(ntpStdHep->mc==mctruth->index && ntpStdHep->IstHEP==3 &&
00872          !(abs(ntpStdHep->IdHEP)==15)){  //not a tau lepton
00873         hfs = ntpStdHep->IdHEP;
00874         break;
00875       }
00876     }
00877     hfs = hfs%1000;
00878   }
00879   else {
00880     for(int i=0;i<nStdHep;i++){
00881       NtpMCStdHep *ntpStdHep = dynamic_cast<NtpMCStdHep *>(stdhepArray[i]);
00882       if(ntpStdHep->mc==mctruth->index && ntpStdHep->IstHEP==3){
00883         hfs = ntpStdHep->IdHEP;
00884         break;
00885       }
00886     }
00887     hfs = hfs%1000;
00888   }
00889   return hfs;
00890 }

int MNtpModule::Initial_state ( NtpMCTruth mctruth,
TClonesArray &  stdhepArray 
) [private]

Definition at line 463 of file MNtpModule.cxx.

References NtpMCStdHep::IdHEP, NtpMCTruth::index, NtpMCStdHep::IstHEP, and NtpMCStdHep::mc.

Referenced by Ana().

00464                                                         {
00465   //copied from MadQuantities::Initial_state
00466   Int_t initial_state=0;
00467   Int_t nStdHep = stdhepArray.GetEntries();
00468 
00469   int protneut = -1;  // 0 = neutron, 1 = proton, 2 = N, 3 = A
00470   int nubarnu = 0;    // +1 = neutrino, -1 = antineutrino
00471 
00472   for(int i=0;i<nStdHep;i++){
00473  
00474     NtpMCStdHep *ntpStdHep = dynamic_cast<NtpMCStdHep *>(stdhepArray[i]);
00475 
00476     if(ntpStdHep->mc==mctruth->index){
00477 
00478       if(ntpStdHep->IstHEP==0){  //initial state particle
00479         if(abs(ntpStdHep->IdHEP)==12 ||
00480            abs(ntpStdHep->IdHEP)==14 ||
00481            abs(ntpStdHep->IdHEP)==16){   //(anti)neutrino
00482           nubarnu = ntpStdHep->IdHEP/abs(ntpStdHep->IdHEP);  //get sign
00483         }
00484       }
00485       if(ntpStdHep->IstHEP==11){    //target nucleon
00486         if(ntpStdHep->IdHEP==2212) protneut = 1;  //proton
00487         else if(ntpStdHep->IdHEP==2112) protneut = 0;  //neutron
00488         else if(abs(ntpStdHep->IdHEP)>1000000000) protneut = 2;  //nucleus
00489         else protneut = 3; //atom - probably never get here since IdHEP A>N?
00490       }
00491     }
00492   }
00493 
00494   if(protneut==1 && nubarnu==1)  initial_state=1;  //p + v
00495   if(protneut==0 && nubarnu==1)  initial_state=2;  //n + v
00496   if(protneut==1 && nubarnu==-1) initial_state=3;  //p + vbar
00497   if(protneut==0 && nubarnu==-1) initial_state=4;  //n + vbar
00498   if(protneut==2 && nubarnu==1)  initial_state=5;  //N + v
00499   if(protneut==3 && nubarnu==1)  initial_state=6;  //A + v
00500   if(protneut==2 && nubarnu==-1) initial_state=7;  //N + vbar
00501   if(protneut==3 && nubarnu==-1) initial_state=8;  //A + vbar
00502 
00503   return initial_state;
00504 }

bool MNtpModule::LoadLargestShowerFromEvent ( NtpSREvent event,
NtpStRecord record,
int &  shwidx 
) [private]

Definition at line 1009 of file MNtpModule.cxx.

References fMNtp, SntpHelpers::GetShower(), NtpSRStripPulseHeight::gev, NtpSRShower::ph, NtpSREvent::shw, and MNtp::shwtotgev.

Referenced by Ana().

01012 {
01013   bool status = false;
01014   float largest_ph = 0;
01015   for (int i=0; i<event->nshower; i++){
01016     NtpSRShower *shower = SntpHelpers::GetShower(event->shw[i],record);
01017     if (!shower) continue;
01018     if (shower->ph.gev>largest_ph){
01019       shwidx = event->shw[i];
01020       largest_ph = shower->ph.gev;
01021     }
01022     if (fMNtp->shwtotgev<=0) fMNtp->shwtotgev = shower->ph.gev;
01023     else fMNtp->shwtotgev += shower->ph.gev;
01024   }
01025   if (largest_ph>0){
01026     status = true;
01027   }
01028   return status;
01029 }

bool MNtpModule::LoadLargestTrackFromEvent ( NtpSREvent event,
NtpStRecord record,
int &  trkidx 
) [private]

Definition at line 942 of file MNtpModule.cxx.

References NtpSRPlane::beg, NtpSRPlane::end, SntpHelpers::GetTrack(), NtpSRTrack::plane, and NtpSREvent::trk.

Referenced by Ana().

00945 {
00946   bool status = false;
00947   float longest_trk = 0;
00948   for (int i=0; i<event->ntrack; i++){
00949     NtpSRTrack *track = SntpHelpers::GetTrack(event->trk[i],record);
00950     if (!track) continue;
00951     int trklen = abs(track->plane.end-track->plane.beg)+1;
00952     if (trklen>longest_trk){
00953       trkidx = event->trk[i];
00954       longest_trk = trklen;
00955     }
00956   }
00957   if (longest_trk>0){
00958     status = true;
00959   }
00960   return status;
00961 }

bool MNtpModule::LoadShowerAtTrackVertex ( NtpSREvent event,
NtpStRecord record,
int  trkidx,
int &  shwidx 
) [private]

Definition at line 965 of file MNtpModule.cxx.

References fMNtp, SntpHelpers::GetShower(), SntpHelpers::GetTrack(), NtpSRStripPulseHeight::gev, NtpSRShowerPulseHeight::linCCgev, NtpSRShower::ph, NtpSREvent::shw, NtpSRShower::shwph, MNtp::shwtotgev, NtpSRShower::vtx, NtpSRTrack::vtx, and NtpSRVertex::z.

Referenced by Ana().

00969 {
00970   bool status = false;
00971   NtpSRTrack *track = SntpHelpers::GetTrack(trkidx,record);
00972   if (!track) return false;
00973   float closest_ph = 0.;
00974   float closest_vtx = 1000.;
00975   const float close_enough = 7*0.06; //about 1 interaction length
00976   for (int i=0; i<event->nshower; i++){
00977     NtpSRShower *shower = SntpHelpers::GetShower(event->shw[i],record);
00978     if (!shower) continue;
00979     if (fMNtp->shwtotgev<=0) fMNtp->shwtotgev = shower->ph.gev;
00980     else fMNtp->shwtotgev += shower->ph.gev;
00981     float vtx_sep = TMath::Abs(shower->vtx.z-track->vtx.z);
00982 
00983     //if this is the closest shower
00984     if (vtx_sep<closest_vtx){
00985       shwidx = event->shw[i];
00986       closest_vtx = vtx_sep;
00987       if (shower->shwph.linCCgev>0) closest_ph = shower->shwph.linCCgev;
00988       else closest_ph = shower->ph.gev;
00989     }
00990     // if this isn't the closest shower but it's within close_enough m of the
00991     // track vertex and has a larger pulseheight
00992     else if(vtx_sep<close_enough && shower->ph.gev>closest_ph){
00993       shwidx = event->shw[i];
00994       closest_vtx = vtx_sep;
00995 
00996       if (shower->shwph.linCCgev>0) closest_ph = shower->shwph.linCCgev;
00997       else closest_ph = shower->ph.gev;
00998     }
00999   }
01000   if (closest_vtx<close_enough){
01001     status = true;
01002   }
01003   else shwidx = -1;
01004   return status;
01005 }

int MNtpModule::Nucleus ( NtpMCTruth mctruth  )  [private]

Definition at line 508 of file MNtpModule.cxx.

References MuELoss::a, NtpMCTruth::a, and NtpMCTruth::z.

Referenced by Ana().

00508                                           {
00509 
00510   Int_t z=int(mctruth->z);
00511   Int_t a=int(mctruth->a);
00512   Int_t nucleus=1;
00513 
00514   switch (z) {
00515     //  case 1:
00516     //nucleus=0;   // free nucleon
00517     //break;
00518   case 1:
00519     switch (a) {
00520     case 1:
00521       nucleus=256;   // hydrogen1
00522       break;
00523     case 2:
00524       nucleus=257;   // hydrogen2
00525       break;
00526     case 3:
00527       nucleus=258;   // hydrogen2
00528       break;
00529     default:
00530       nucleus=256;   // hydrogen1
00531       break;
00532     }
00533     break;
00534   case 6:
00535     switch (a) {
00536     case 11:
00537       nucleus=274; // carbon11   
00538       break;
00539     case 12:
00540       nucleus=275; // carbon12
00541       break;
00542     case 13:
00543       nucleus=276; // carbon13
00544       break;
00545     case 14:
00546       nucleus=277; // carbon14
00547       break;
00548     default:
00549       nucleus=275; // carbon12
00550       break;
00551     }
00552     break;
00553   case 7:
00554     switch (a) {
00555     case 13:
00556       nucleus=278; // nitrogen13   
00557       break;
00558     case 14:
00559       nucleus=279; // nitrogen14
00560       break;
00561     case 15:
00562       nucleus=280; // nitrogen15
00563       break;
00564     case 16:
00565       nucleus=281; // nitrogen16
00566       break;
00567     case 17:
00568       nucleus=282; // nitrogen17
00569       break;
00570     default:
00571       nucleus=279; // nitrogen14
00572       break;
00573     }
00574     break;
00575   case 8:
00576     switch (a) {
00577     case 15:
00578       nucleus=283;   // oxygen15
00579       break;
00580     case 16:
00581       nucleus=284;   // oxygen16
00582       break;
00583     case 17:
00584       nucleus=285;   // oxygen17
00585       break;
00586     case 18:
00587       nucleus=286;   // oxygen18
00588       break;
00589     default:
00590       nucleus=284;   // oxygen16
00591       break;
00592     }
00593     break;
00594   case 13:
00595     switch (a) {
00596     case 26:
00597       nucleus=303;   // aluminium26
00598       break;
00599     case 27:
00600       nucleus=304;   // aluminium27
00601       break;
00602     case 28:
00603       nucleus=305;   // aluminium28
00604       break;
00605     case 29:
00606       nucleus=306;   // aluminium29
00607       break;
00608     default:
00609       nucleus=304;   // aluminium27
00610       break;
00611     }
00612     break;
00613   case 14:
00614     switch (a) {
00615     case 27:
00616       nucleus=307;   // silicon27
00617       break;
00618     case 28:
00619       nucleus=308;   // silicon28
00620       break;
00621     case 29:
00622       nucleus=309;   // silicon29
00623       break;
00624     case 30:
00625       nucleus=310;   // silicon30
00626       break;
00627     default:         
00628       nucleus=308;   // silicon28
00629       break;
00630     }
00631     break;
00632   case 15:
00633     switch (a) {
00634     case 30:
00635       nucleus=311;   //phosphorus30
00636       break;
00637     case 31:
00638       nucleus=312;   //phosphorus31
00639       break;
00640     case 32:
00641       nucleus=313;   //phosphorus32
00642       break;
00643     case 33:
00644       nucleus=314;   //phosphorus33
00645       break;
00646     default:
00647       nucleus=312;   //phosphorus31
00648       break;
00649     }
00650     break;
00651   case 16:
00652     switch (a) {
00653     case 31:
00654       nucleus=315;   //sulphur31
00655       break;
00656     case 32:
00657       nucleus=316;   //sulphur32
00658       break;
00659     case 33:
00660       nucleus=317;   //sulphur33
00661       break;
00662     case 34:
00663       nucleus=318;   //sulphur34
00664       break;
00665     case 35:
00666       nucleus=319;   //sulphur35
00667       break;
00668     case 36:
00669       nucleus=320;   //sulphur36
00670       break;
00671     default:
00672       nucleus=316;   //sulphur32
00673       break;
00674     }
00675     break;
00676   case 22:
00677     switch (a) {
00678     case 45:
00679       nucleus=347;   //titanium45
00680       break;
00681     case 46:
00682       nucleus=348;   //titanium46
00683       break;
00684     case 47:
00685       nucleus=349;   //titanium47
00686       break;
00687     case 48:
00688       nucleus=350;   //titanium48
00689       break;
00690     case 49:
00691       nucleus=351;   //titanium49
00692       break;
00693     case 50:
00694       nucleus=352;   //titanium50
00695       break;
00696     default:
00697       nucleus=350;   //titanium48
00698       break;
00699     }
00700     break;
00701   case 23:
00702     switch (a) {
00703     case 49:
00704       nucleus=353;   //vanadium49
00705       break;
00706     case 50:
00707       nucleus=354;   //vanadium50
00708       break;
00709     case 51:
00710       nucleus=355;   //vanadium51
00711       break;
00712     case 52:
00713       nucleus=356;   //vanadium52
00714       break;
00715     case 53:
00716       nucleus=357;   //vanadium53
00717       break;
00718     default:
00719       nucleus=355;   //vanadium51
00720       break;
00721     }
00722     break;
00723   case 24:
00724     switch (a) {
00725     case 49:
00726       nucleus=358;   //chromium49
00727       break;
00728     case 50:
00729       nucleus=359;   //chromium50
00730       break;
00731     case 51:
00732       nucleus=360;   //chromium51
00733       break;
00734     case 52:
00735       nucleus=361;   //chromium52
00736       break;
00737     case 53:
00738       nucleus=362;   //chromium53
00739       break;
00740     case 54:
00741       nucleus=363;   //chromium54
00742       break;
00743     default:
00744       nucleus=361;   //chromium52
00745       break;
00746     }
00747     break;
00748   case 25:
00749     switch (a) {
00750     case 53:
00751       nucleus=364;   //manganese53
00752       break;
00753     case 54:
00754       nucleus=365;   //manganese54
00755       break;    
00756     case 55:
00757       nucleus=366;   //manganese55
00758       break;    
00759     case 56:
00760       nucleus=367;   //manganese56
00761       break;    
00762     case 57:
00763       nucleus=368;   //manganese57
00764       break;    
00765     default:
00766       nucleus=366;   //manganese55
00767       break;
00768     }
00769     break;
00770   case 26:
00771     switch (a) {
00772     case 53:
00773       nucleus=369;   //iron53
00774       break;
00775     case 54:
00776       nucleus=370;   //iron54
00777       break;
00778     case 55:
00779       nucleus=371;   //iron55
00780       break;
00781     case 56:
00782       nucleus=372;   //iron56
00783       break;
00784     case 57:
00785       nucleus=373;   //iron57
00786       break;
00787     case 58:
00788       nucleus=374;   //iron58
00789       break;
00790     default:
00791       nucleus=372;   //iron56
00792       break;
00793     }
00794     break;
00795   case 28:
00796     switch (a) {
00797     case 57:
00798       nucleus=382;   //nickel57
00799       break;
00800     case 58:
00801       nucleus=383;   //nickel58
00802       break;
00803     case 59:
00804       nucleus=384;   //nickel59
00805       break;
00806     case 60:
00807       nucleus=385;   //nickel60
00808       break;
00809     case 61:
00810       nucleus=386;   //nickel61
00811       break;
00812     case 62:
00813       nucleus=387;   //nickel62
00814       break;
00815     case 63:
00816       nucleus=388;   //nickel63
00817       break;
00818     case 64:
00819       nucleus=389;   //nickel64
00820       break;
00821     default:
00822       nucleus=383;   //nickel58
00823       break;
00824     }
00825     break;
00826   case 29:
00827     switch (a) {
00828     case 62:
00829       nucleus=390;   //copper62
00830       break;
00831     case 63:
00832       nucleus=391;   //copper63
00833       break;
00834     case 64:
00835       nucleus=392;   //copper64
00836       break;
00837     case 65:
00838       nucleus=393;   //copper65
00839       break;
00840     case 66:
00841       nucleus=394;   //copper66
00842       break;
00843     case 67:
00844       nucleus=395;   //copper67
00845       break;
00846     default:
00847       nucleus=392;   //copper64
00848       break;
00849     }
00850     break;
00851   default:
00852     nucleus=1;   // unknown
00853     break;
00854   }
00855 
00856   return nucleus;
00857 }

bool MNtpModule::PlanTriggerSim ( NtpSREvent event,
NtpStRecord record 
) [private]

Definition at line 894 of file MNtpModule.cxx.

References fM, fN, VHS::GetStrip(), NtpSRStrip::ndigit, NtpSRPulseHeight::pe, NtpSRStrip::ph1, NtpSRStrip::plane, and NtpSREvent::stp.

Referenced by Ana().

00896 {
00897   int nstp = event->nstrip; //no. of strips
00898 
00899   std::vector<Int_t> planeMap(600);
00900   for(UInt_t i = 0;i<planeMap.size();i++) planeMap[i]=0;
00901 
00902   Int_t totplanes = 0;
00903   Int_t firstPlane = 9999;
00904   Int_t lastPlane = -9999;
00905 
00906   for (int istp = 0; istp<nstp; istp++){
00907     NtpSRStrip *strip = SntpHelpers::GetStrip(event->stp[istp],record);
00908     Int_t plane = strip->plane ;
00909 
00910     if(plane > (Int_t)planeMap.size()) continue;
00911     if(strip->ph1.pe == 0 || strip->ndigit==0) continue ;
00912     
00913     // Add plane to map.
00914     if(planeMap[plane] ==0){
00915       totplanes++;
00916       if(plane<firstPlane) firstPlane = plane;
00917       if(plane>lastPlane) lastPlane = plane;
00918     }
00919     
00920     planeMap[plane]++;
00921   }
00922 
00923   if(totplanes < fN) return false;
00924   
00925   for(int i = firstPlane; i <= lastPlane; i++) {
00926     Int_t nPlane = 0;
00927     
00928     for (Int_t j=0; (j<fM) && (i+j<=lastPlane); j++) {
00929       if ( planeMap[i+j]>0 ) nPlane++;
00930     }
00931     if (nPlane>=fN) {
00932       return true;
00933     }
00934     if (i+fM>lastPlane) break; // We're done... stop looking.
00935   }
00936   
00937   return false;
00938 }


Member Data Documentation

Definition at line 68 of file MNtpModule.h.

Referenced by Ana().

Definition at line 65 of file MNtpModule.h.

Referenced by Ana().

Definition at line 64 of file MNtpModule.h.

Referenced by Ana().

Definition at line 78 of file MNtpModule.h.

Referenced by Ana().

std::string MNtpModule::fFile [private]

Definition at line 56 of file MNtpModule.h.

Referenced by BeginJob(), and Config().

Int_t MNtpModule::fM [private]

Definition at line 58 of file MNtpModule.h.

Referenced by Ana(), and PlanTriggerSim().

MNtp* MNtpModule::fMNtp [private]
Int_t MNtpModule::fN [private]

Definition at line 57 of file MNtpModule.h.

Referenced by Ana(), and PlanTriggerSim().

TFile* MNtpModule::fNtpFile [private]

Definition at line 62 of file MNtpModule.h.

Referenced by BeginJob(), and EndJob().

TTree* MNtpModule::fNtuple [private]

Definition at line 63 of file MNtpModule.h.

Referenced by Ana(), BeginJob(), and EndJob().

Definition at line 79 of file MNtpModule.h.

Referenced by Ana().

bool MNtpModule::IsTriggerOn [private]

Definition at line 59 of file MNtpModule.h.

Referenced by Ana(), and Config().

int MNtpModule::numfiles [private]

Definition at line 75 of file MNtpModule.h.

Referenced by BeginFile().

Definition at line 72 of file MNtpModule.h.

Referenced by FillTrkInfo(), and MNtpModule().

Definition at line 71 of file MNtpModule.h.

Referenced by FillTrkInfo(), and MNtpModule().


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

Generated on 2 Nov 2017 for loon by  doxygen 1.6.1