MCAnalysis Class Reference

#include <MCAnalysis.h>

Inheritance diagram for MCAnalysis:
JobCModule

List of all members.

Public Member Functions

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

Private Member Functions

void FillTrkInfo (const CandTrackHandle *track)
void FillPlnInfo (const CandDeMuxDigitListHandle *digitlist)
int GetMaxTrkPlane ()

Private Attributes

std::string fNameListIn
Int_t fMakePictures
TFile * fFileOut
TTree * fTreeOut
Int_t fRun
Int_t fSnarl
Int_t fNu
Int_t fAction
Float_t fENu
Float_t fEMu
Int_t fNMuonDig
Int_t fNMuonDigRetained
Int_t fNShwDig
Int_t fNShwDigRetained
Int_t fNShwDigAtVtx
Int_t fNShwDigRetainedAtVtx
Float_t fNShwPE
Float_t fNShwPERetained
Float_t fNShwPEAtVtx
Float_t fNShwPERetainedAtVtx
Int_t fNRetained
Int_t fNRetainedMuon
Int_t fNRetainedShw
Int_t fNRetainedBoth
Float_t fPERetained
Float_t fPERetainedMuon
Float_t fPERetainedShw
Float_t fPERetainedBoth
Int_t fNRejected
Int_t fNRejectedMuon
Int_t fNRejectedShw
Int_t fNRejectedBoth
Int_t fRShwN
Float_t fRShwE [500]
Int_t fRShwDist [500]
Int_t fRShwMuon [500]
Int_t fRShwTrk [500]
Int_t fRShwXTalk [500]
Int_t fRShwHitType [500]
Int_t fRShwHitStatus [500]
Int_t fRShwHitMother [500]
Int_t fTrkMIPCalibOkay
TH1F * hMIPTrueMuon
TH1F * hMIPMissTrack
TH1F * hMIPMixed
TH1F * hMIPTrueMuon2
TH1F * hMIPMissTrack2
TH1F * hMIPMixed2
TCanvas * fCanvas
TH2F * fUZAxis
TH2F * fVZAxis
vector< TMarker > fUHits
vector< TMarker > fVHits
TCanvas * fCanvas2
TH1F * fUZTrkQ
TH1F * fVZTrkQ
Float_t fTrkQ [500]
Float_t fTrkMIPDcosz [500]
Float_t fPlnQ [500]
set< int > fTrkStrips

Detailed Description

Definition at line 26 of file MCAnalysis.h.


Constructor & Destructor Documentation

MCAnalysis::MCAnalysis (  ) 

(Document me!)

Definition at line 43 of file MCAnalysis.cxx.

00043                        : 
00044   fNameListIn(""),
00045 
00046   fMakePictures(0)
00047 {
00051 }

MCAnalysis::~MCAnalysis (  ) 

(Document me!)

Definition at line 54 of file MCAnalysis.cxx.

00055 {
00059 }


Member Function Documentation

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

Implement this for read only access to the MomNavigator

Reimplemented from JobCModule.

Definition at line 168 of file MCAnalysis.cxx.

References digit(), Draw(), fAction, fCanvas, fCanvas2, fEMu, fENu, FillPlnInfo(), FillTrkInfo(), CandRecord::FindCandHandle(), RecDataRecord< T >::FindComponent(), fMakePictures, fNMuonDig, fNMuonDigRetained, fNRejected, fNRejectedBoth, fNRejectedMuon, fNRejectedShw, fNRetained, fNRetainedBoth, fNRetainedMuon, fNRetainedShw, fNShwDig, fNShwDigAtVtx, fNShwDigRetained, fNShwDigRetainedAtVtx, fNShwPE, fNShwPEAtVtx, fNShwPERetained, fNShwPERetainedAtVtx, fNu, fPERetained, fPERetainedBoth, fPERetainedMuon, fPERetainedShw, fRShwDist, fRShwE, fRShwHitMother, fRShwHitStatus, fRShwHitType, fRShwMuon, fRShwN, fRShwTrk, fRShwXTalk, fRun, fSnarl, fTreeOut, fTrkMIPCalibOkay, fTrkMIPDcosz, fTrkStrips, fUHits, fUZAxis, fUZTrkQ, fVHits, fVZAxis, fVZTrkQ, Truthifier::GetBiggestHit(), CandHandle::GetDaughterIterator(), CandDeMuxDigitHandle::GetDeMuxDigitFlagWord(), MomNavigator::GetFragment(), RecMinos::GetHeader(), DigiSignal::GetHit(), GetMaxTrkPlane(), CandHandle::GetNDaughters(), DigiSignal::GetNumberOfHits(), CandHeader::GetRun(), Truthifier::GetSignal(), VldContext::GetSimFlag(), CandHeader::GetSnarl(), DigiSignal::GetTruth(), RecMinos::GetVldContext(), header, hMIPMissTrack, hMIPMissTrack2, hMIPMixed, hMIPMixed2, hMIPTrueMuon, hMIPTrueMuon2, REROOT_NeuKin::IAction(), REROOT_NeuKin::INu(), TruthHelper::IsMuon(), TruthHelper::IsNeutrino(), Truthifier::IsValid(), Msg::kError, JobCResult::kFailed, DigiSignal::kGenuine, SimFlag::kMC, JobCResult::kPassed, CalDigitType::kPE, PlaneView::kU, Msg::kWarning, CandDeMuxDigit::kXTalk, MSG, REROOT_NeuKin::P4Mu1(), REROOT_NeuKin::P4Neu(), and DigiScintHit::TrackId().

00169 {  
00170   //
00171   //Set up all the helper objects we need
00172   //
00173   CandRecord* record = dynamic_cast<CandRecord*>(mom->GetFragment("CandRecord"));
00174   if(!record){
00175     MSG("MCAna",Msg::kError) << " Unable to find CandRecord in mom !!! " << endl;
00176     return JobCResult::kFailed;
00177   }
00178   const VldContext &vldc = *(record->GetVldContext());      
00179   if(vldc.GetSimFlag()!=SimFlag::kMC){
00180     MSG("MCAna",Msg::kWarning) << " Event not a MC event"<<endl;
00181     return JobCResult::kPassed;
00182   }
00183   Truthifier* truth = new Truthifier(mom);
00184   TruthHelper* beer = new TruthHelper(mom); 
00185   if(!(truth->IsValid())){
00186     MSG("MCAna",Msg::kError) << " The THRUTH is not valid! "<<endl;
00187     return JobCResult::kFailed;
00188   }
00189   //
00190   //Get pointers to all the reco object we need
00191   //
00192   const CandEventListHandle * eventlist = dynamic_cast<const CandEventListHandle*>(record->FindCandHandle("CandEventListHandle"));
00193   if(!eventlist || eventlist->GetNDaughters()!=1){
00194     MSG("MCAna",Msg::kWarning) << " Rejecting event as it has no events " << endl;
00195     return JobCResult::kFailed;
00196   }
00197 
00198 
00199   TIter event_iter(eventlist->GetDaughterIterator());
00200   const CandEventHandle* event = dynamic_cast<const CandEventHandle*>(event_iter());
00201   assert(event);
00202 
00203   //build a list of digits that are in the track
00204   const CandTrackHandle* track = dynamic_cast<const CandTrackHandle*>(event->GetPrimaryTrack());
00205   assert(track);
00206   
00207   const CandDigitListHandle* rawdigitlist = dynamic_cast<const CandDigitListHandle*>(record->FindCandHandle("CandDigitListHandle", "altdemux"));
00208   const CandDigitListHandle* strippeddigitlist = dynamic_cast<const CandDigitListHandle*>(record->FindCandHandle("CandDigitListHandle", "stripdigitlist"));
00209   
00210   if(rawdigitlist==NULL || strippeddigitlist==NULL){
00211     MSG("MCAna",Msg::kError) << " Unable to find Merged digit list (" << strippeddigitlist<<")  or RawDigitList (" << rawdigitlist <<")"<< endl;
00212     return JobCResult::kFailed;
00213   }
00214   
00215   const CandHeader* header  = dynamic_cast<const CandHeader*>(record->GetHeader());
00216   //
00217   //Clear the display objects
00218   //
00219   fUHits.clear();
00220   fVHits.clear();
00221   fUZTrkQ->Reset();
00222   fVZTrkQ->Reset();
00223   
00224   //
00225   //Only perform this analysis on mucc
00226   //
00227   SimSnarlRecord *simsnarl = dynamic_cast<SimSnarlRecord*>(mom->GetFragment("SimSnarlRecord"));
00228   if (!simsnarl){
00229     MSG("MCAna", Msg::kError)<< " Faied to get the truth (SimSnarl) " <<endl; 
00230     return JobCResult::kFailed;
00231   }
00232   const TClonesArray* ctca = dynamic_cast<const TClonesArray*>(simsnarl->FindComponent("TClonesArray","StdHep"));
00233   const TClonesArray* neukinarray = dynamic_cast<const TClonesArray*>(simsnarl->FindComponent("TClonesArray","NeuKinList"));
00234   const REROOT_NeuKin* neukin = dynamic_cast<const REROOT_NeuKin*>(neukinarray->At(0));
00235   
00236   assert(neukin);
00237   const bool numucc = (abs(neukin->INu())==14 && abs(neukin->IAction()!=0) );
00238    
00239   //
00240   //build a list of digits that are retained
00241   //
00242   set<Int_t> retainedDigitIds;
00243   TIter strippeddigititer(strippeddigitlist->GetDaughterIterator());
00244   while (CandDigitHandle *digit = dynamic_cast<CandDigitHandle*>(strippeddigititer())) { 
00245     if(retainedDigitIds.count(digit->GetRawDigitIndex())!=0){
00246       MSG("MCAna", Msg::kError)<< " Mulitple RAWDIGITID entries " <<endl; 
00247     }else{
00248       if(digit->GetQieErrorBits()==0){
00249         retainedDigitIds.insert(digit->GetRawDigitIndex());
00250       }else{
00251         retainedDigitIds.insert(-1*digit->GetRawDigitIndex());
00252       }
00253     }
00254   }
00255 
00256   //
00257   //For each digit decide if it 'should' have been removed
00258   //
00259   fRun   = header->GetRun();
00260   fSnarl = header->GetSnarl();
00261   fNu    =  neukin->INu();
00262   fAction = neukin->IAction();
00263   fENu   = fabs((neukin->P4Neu())[3]);
00264   fEMu   = fabs((neukin->P4Mu1())[3]);
00265   
00266   fNMuonDig         = 0;
00267   fNMuonDigRetained = 0; 
00268 
00269   fNShwDig          = 0;
00270   fNShwDigRetained  = 0;
00271 
00272   fNShwDigAtVtx          = 0;
00273   fNShwDigRetainedAtVtx  = 0;
00274 
00275   fNShwPE = 0 ;
00276   fNShwPERetained =0;  
00277 
00278   fNShwPEAtVtx = 0 ;
00279   fNShwPERetainedAtVtx =0;  
00280   
00281   fNRetained  = 0 ;
00282   fNRetainedMuon = 0;
00283   fNRetainedShw = 0;
00284   fNRetainedBoth = 0;
00285 
00286   fPERetained  = 0 ;
00287   fPERetainedMuon = 0;
00288   fPERetainedShw = 0;
00289   fPERetainedBoth = 0;
00290 
00291   fNRejected = 0;
00292   fNRejectedMuon = 0;
00293   fNRejectedShw = 0;
00294   fNRejectedBoth = 0;
00295 
00296   fRShwN = 0 ;
00297   for(int i=0; i<500; i++){
00298     fRShwE[i] = 0;
00299     fRShwDist[i] = 0;
00300     fRShwMuon[i] = 0;
00301     fRShwTrk[i] = 0;    
00302     fRShwHitType[i] = 0;    
00303     fRShwHitStatus[i] = 0;    
00304     fRShwHitMother[i] = 0;    
00305   }
00306   
00307   //
00308   //Local variables 
00309   //
00310   int lNRejShw = 0;
00311   int lNRejShwMaxTrk = 0 ;
00312   int lNRejShwFakeTrk = 0;
00313   int lNRejShwMix = 0 ;
00314 
00315   //
00316   //build local arrays
00317   //
00318   FillPlnInfo((CandDeMuxDigitListHandle*)rawdigitlist);
00319   FillTrkInfo(track);
00320   
00321   const int maxtrkplane = GetMaxTrkPlane();
00322 
00323   int minplane = 500;
00324   int maxplane = 0;
00325   
00326 
00327   bool plane_track[500];
00328   bool plane_muon[500];
00329   bool plane_shw[500];
00330   for(int i =0; i<500; i++){
00331     plane_track[i] = plane_muon[i] = plane_shw[i] = false;
00332   }
00333   
00334   TIter rawdigititer(rawdigitlist->GetDaughterIterator());
00335   
00336   while (CandDigitHandle *digit = dynamic_cast<CandDigitHandle*>(rawdigititer())) {
00337 
00338     const int planeno = digit->GetPlexSEIdAltL().GetBestSEId().GetPlane();
00339     const int stripno = digit->GetPlexSEIdAltL().GetBestSEId().GetStrip();    
00340     const CandDigitHandle tcdh = *digit;      
00341     const CandDeMuxDigitHandle* demuxdigit = dynamic_cast<const CandDeMuxDigitHandle*>(digit);
00342     assert(demuxdigit);
00343     bool  is_shower = 0;
00344     bool  is_muon = 0;
00345     bool is_physics = 0;
00346 
00347     const bool is_retained = (retainedDigitIds.count(digit->GetRawDigitIndex())>0);        
00348     bool is_track = (fTrkStrips.count(planeno*192+stripno)>0);
00349     const bool is_recoed_xtalk = (demuxdigit->GetDeMuxDigitFlagWord()&CandDeMuxDigit::kXTalk); 
00350     const bool is_retained_scaled = (retainedDigitIds.count(-1*digit->GetRawDigitIndex())>0);        
00351 
00352     const DigiSignal * signal = truth->GetSignal(tcdh);  
00353 
00354     const DigiScintHit* biggest_hit = truth->GetBiggestHit(tcdh);
00355     TParticle* biggest_part  = NULL;
00356     if(biggest_hit!=NULL){
00357       Int_t biggest_track = abs(biggest_hit->TrackId());
00358       biggest_part = dynamic_cast<TParticle*>((*ctca)[biggest_track]);      
00359     }
00360 
00361     if(signal!=NULL){      
00362       is_physics = ((signal->GetTruth() & DigiSignal::kGenuine)!=0);      
00363       for (UInt_t i=0;i<signal->GetNumberOfHits();i++){
00364         const DigiScintHit * hit = signal->GetHit(i);
00365         if(hit){
00366           Int_t track = abs(hit->TrackId());
00367           TParticle* part = dynamic_cast<TParticle*>((*ctca)[track]);
00368           if(!beer->IsNeutrino(part) && !beer->IsMuon(part)){
00369             is_shower=true;         
00370           }
00371           if(!beer->IsNeutrino(part) && beer->IsMuon(part)){
00372             is_muon=true;           
00373           }     
00374         }
00375       }         
00376     }
00377     
00378     
00379     if(is_physics){
00380       if(is_shower){
00381         fNShwDig++;
00382         fNShwPE+=digit->GetCharge(CalDigitType::kPE);   
00383         if(planeno<maxtrkplane){
00384           fNShwDigAtVtx++;
00385           fNShwPEAtVtx+=digit->GetCharge(CalDigitType::kPE);    
00386         }
00387         if(is_retained || is_retained_scaled){
00388           fNShwDigRetained++;
00389           fNShwPERetained+=digit->GetCharge(CalDigitType::kPE); 
00390           if(planeno<maxtrkplane){
00391             fNShwDigRetainedAtVtx++;
00392             fNShwPERetainedAtVtx+=digit->GetCharge(CalDigitType::kPE);  
00393           }
00394         }else{
00395           if(fRShwN<1000){
00396             fRShwE[fRShwN]=digit->GetCharge(CalDigitType::kPE);
00397             fRShwDist[fRShwN]=maxtrkplane - planeno;
00398             fRShwMuon[fRShwN]=(int)is_muon;
00399             fRShwTrk[fRShwN]=(int)(is_track && !is_muon);
00400             fRShwXTalk[fRShwN]= (int)(is_recoed_xtalk);
00401             if(biggest_part!=NULL){
00402               fRShwHitType[fRShwN] = biggest_part->GetPdgCode();            
00403               fRShwHitStatus[fRShwN] = biggest_part->GetStatusCode();
00404               fRShwHitMother[fRShwN] = biggest_part->GetMother(0);          
00405             }
00406             fRShwN++;
00407           }
00408           lNRejShw++;
00409           if(planeno>maxtrkplane) lNRejShwMaxTrk++;
00410           if(is_muon) lNRejShwMix++;
00411           if(!is_muon && is_track) lNRejShwFakeTrk++;
00412         }
00413       }//is shower
00414       if(is_muon){
00415         fNMuonDig++;
00416         if(is_retained || is_retained_scaled) fNMuonDigRetained++;
00417       }//is muon
00418       if(is_retained || is_retained_scaled){
00419         fNRetained++;
00420         fPERetained+=digit->GetCharge(CalDigitType::kPE);       
00421         if(is_muon){
00422           fNRetainedMuon++; 
00423           fPERetainedMuon+=digit->GetCharge(CalDigitType::kPE); 
00424         }
00425         if(is_shower){
00426           fNRetainedShw++;
00427           fPERetainedShw+=digit->GetCharge(CalDigitType::kPE);  
00428         }
00429         if(is_muon && is_shower){
00430           fNRetainedBoth++;
00431           fPERetainedBoth+=digit->GetCharge(CalDigitType::kPE); 
00432         }
00433       }else{
00434         fNRejected++;
00435         if(is_muon && !is_shower)fNRejectedMuon++;
00436         if(!is_muon && is_shower) fNRejectedShw++;
00437         if(is_muon && is_shower) fNRejectedBoth++;
00438       }     
00439       //
00440       //Decide what plane type this is
00441       //
00442       if(is_track){
00443         plane_track[planeno] = true;
00444         if(is_muon) plane_muon[planeno] = true;
00445         if(is_shower) plane_shw[planeno] = true;
00446       }
00447 
00448       //                                        
00449       //Drawing
00450       //
00451       TMarker marker(planeno, stripno,25);    
00452       if( is_muon && !is_shower) marker.SetMarkerStyle(28);
00453       if(!is_muon &&  is_shower) marker.SetMarkerStyle(24);
00454       if( is_muon &&  is_shower) marker.SetMarkerStyle(26);
00455       if(is_retained)marker.SetMarkerColor(1);
00456       else if(is_retained_scaled)marker.SetMarkerColor(4);
00457       else marker.SetMarkerColor(2);    
00458       if(digit->GetPlexSEIdAltL().GetBestSEId().GetPlaneView()==PlaneView::kU){
00459         fUHits.push_back(marker);
00460       }else{
00461         fVHits.push_back(marker);
00462       }
00463       if(is_physics && is_shower){
00464         if(planeno>maxplane) maxplane = planeno;
00465         if(planeno<minplane) minplane = planeno;      
00466       }           
00467     }//is_physics
00468 
00469   }
00470   //
00471   //fill info about the event
00472   //
00473   int ntrkcalibokay = 0;
00474   for(int i=0; i<500; i++){
00475     if(plane_track[i]){
00476       if(fTrkMIPDcosz[i]!=0){
00477         ntrkcalibokay++;
00478         if(plane_muon[i] && !plane_shw[i]) hMIPTrueMuon ->Fill(fTrkMIPDcosz[i]);
00479         if(plane_muon[i] &&  plane_shw[i]) hMIPMixed    ->Fill(fTrkMIPDcosz[i]);
00480         if(!plane_muon[i])                 hMIPMissTrack->Fill(fTrkMIPDcosz[i]);
00481         if(i<maxtrkplane){
00482           if(plane_muon[i] && !plane_shw[i]) hMIPTrueMuon2 ->Fill(fTrkMIPDcosz[i]);
00483           if(plane_muon[i] &&  plane_shw[i]) hMIPMixed2    ->Fill(fTrkMIPDcosz[i]);
00484           if(!plane_muon[i])                 hMIPMissTrack2->Fill(fTrkMIPDcosz[i]);     
00485         }
00486       }
00487     }
00488   }
00489   fTrkMIPCalibOkay = (ntrkcalibokay>4);
00490   fTreeOut->Fill();
00491 
00492   if(fMakePictures){
00493     cout << " Event MuCC: " << ((numucc)?"YES":"NO")  <<endl;
00494     cout << " Max Trk Pln: "<< maxtrkplane<<endl;
00495     cout << " Muon Digits   : " << fNMuonDig<<endl;
00496     cout << "    Retained   : " << fNMuonDigRetained / (float)fNMuonDig << " (" << fNMuonDigRetained << ")"<<endl;
00497     cout << "    Rejected   : " << (fNMuonDig - fNMuonDigRetained) / (float)fNMuonDig << " (" << fNMuonDig - fNMuonDigRetained << ")"<<endl;
00498     cout << " Shower Digits : " << fNShwDig <<endl; 
00499     cout << "    Retained   : " << fNShwDigRetained / (float)fNShwDig << " (" << fNShwDigRetained <<")" <<endl;
00500     cout << "    Rejected   : " << (fNShwDig - fNShwDigRetained) / (float)fNShwDig << " (" << fNShwDig - fNShwDigRetained <<")" <<endl;
00501     
00502     
00503     cout << " Rejected    : " << fNRejected<<endl;
00504     cout << "    Muon only: " << fNRejectedMuon/(float)fNRejected << " "  << fNRejectedMuon<<endl;
00505     cout << "    Shw only : " << fNRejectedShw/(float)fNRejected << " "  << fNRejectedShw<<endl;
00506     cout << "    both     : " << fNRejectedBoth/(float)fNRejected << " "  << fNRejectedBoth<<endl;
00507     cout << " Retained    : " << fNRetained<<endl;
00508     cout << "    Muon only: " << fNRetainedMuon/(float)fNRetained << " "  << fNRetainedMuon<<endl;
00509     cout << "    Shw only : " << fNRetainedShw/(float)fNRetained << " "  << fNRetainedShw<<endl;
00510     cout << "    both     : " << fNRetainedBoth/(float)fNRetained << " "  << fNRetainedBoth<<endl;
00511     
00512     cout << " Rejected Shw: " <<endl;
00513     cout << "        Total: " << lNRejShw <<endl;
00514     if(lNRejShw>0){
00515       cout << "       MAXTRK: " << lNRejShwMaxTrk  << "  " <<  lNRejShwMaxTrk / lNRejShw<<endl;
00516       cout << "        MIXED: " << lNRejShwMix << "  " << lNRejShwMix / lNRejShw<<endl;
00517       cout << "      FAKETRK: " << lNRejShwFakeTrk << "  " <<  lNRejShwFakeTrk / lNRejShw<<endl;
00518     }
00519     cout << " Calibrated: " << ((fTrkMIPCalibOkay!=0)? " yes": "NO***") <<endl;
00520     
00521     fCanvas ->Clear();
00522     fCanvas->Divide(2,1);
00523     fCanvas->cd(1);
00524     fUZAxis->GetXaxis()->SetRangeUser(minplane-5, maxplane+10);
00525     fUZAxis->Draw();
00526     for(unsigned int i=0; i<fUHits.size(); i++) fUHits[i].Draw();
00527     fCanvas->cd(2);
00528     fVZAxis->GetXaxis()->SetRangeUser(minplane-5, maxplane+10);
00529     fVZAxis->Draw();
00530     for(unsigned int i=0; i<fVHits.size(); i++) fVHits[i].Draw();
00531     fCanvas->Draw();
00532     fCanvas->Update();
00533     
00534     fCanvas2->Clear();
00535     fCanvas2->Divide(2,1);
00536     fCanvas2->cd(1);
00537     fUZTrkQ->GetXaxis()->SetRangeUser(minplane-5, maxplane+10);
00538     fUZTrkQ->Draw("hist");
00539     fCanvas2->cd(2);
00540     fVZTrkQ->GetXaxis()->SetRangeUser(minplane-5, maxplane+10);
00541     fVZTrkQ->Draw("hist");
00542     fCanvas2->Draw();
00543     fCanvas2->Update();
00544     
00545     char histname[1024];
00546     sprintf(histname, "mc_event_%d_%d.eps", header->GetRun(), header->GetSnarl());
00547     fCanvas->Print(histname);
00548     sprintf(histname, "mc_q_%d_%d.eps", header->GetRun(), header->GetSnarl());
00549     fCanvas2->Print(histname);
00550   }
00551   delete truth;
00552   delete beer;
00553   return JobCResult::kPassed; // kNoDecision, kFailed, etc.
00554 }

void MCAnalysis::BeginJob ( void   )  [virtual]

Implement for notification of begin of job

(Document me!)

Reimplemented from JobCModule.

Definition at line 63 of file MCAnalysis.cxx.

References fAction, fCanvas, fCanvas2, fEMu, fENu, fFileOut, fNMuonDig, fNMuonDigRetained, fNRejected, fNRejectedBoth, fNRejectedMuon, fNRejectedShw, fNRetained, fNRetainedBoth, fNRetainedMuon, fNRetainedShw, fNShwDig, fNShwDigAtVtx, fNShwDigRetained, fNShwDigRetainedAtVtx, fNShwPE, fNShwPEAtVtx, fNShwPERetained, fNShwPERetainedAtVtx, fNu, fPERetained, fPERetainedBoth, fPERetainedMuon, fPERetainedShw, fRShwDist, fRShwE, fRShwHitMother, fRShwHitStatus, fRShwHitType, fRShwMuon, fRShwN, fRShwTrk, fRShwXTalk, fRun, fSnarl, fTreeOut, fTrkMIPCalibOkay, fUZAxis, fUZTrkQ, fVZAxis, fVZTrkQ, hMIPMissTrack, hMIPMissTrack2, hMIPMixed, hMIPMixed2, hMIPTrueMuon, and hMIPTrueMuon2.

00064 {
00068   fCanvas = new TCanvas("mcana", "MC Analysis", 800, 600);
00069   fUZAxis = new TH2F("uzaxis", "", 500, 0, 500, 192, 0, 192);
00070   fVZAxis = new TH2F("vzaxis", "", 500, 0, 500, 192, 0, 192);
00071   fCanvas2 = new TCanvas("mcana2", "MC Analysis (Q)", 800, 600);
00072   fUZTrkQ = new TH1F("trkquz",  "", 500, 0, 500);
00073   fVZTrkQ = new TH1F("trkqvz",  "", 500, 0, 500);
00074   
00075   TDirectory* cdir = gDirectory;
00076   fFileOut = TFile::Open("mc_performance.root", "recreate");
00077   fFileOut->cd();
00078   fTreeOut = new TTree("mcp", "MC performance");
00079   fTreeOut->SetDirectory(fFileOut);
00080   fTreeOut->Branch("run", &fRun, "run/I");
00081   fTreeOut->Branch("snarl", &fSnarl, "snarl/I");
00082   fTreeOut->Branch("nu", &fNu, "nu/I");
00083   fTreeOut->Branch("action", &fAction, "action/I");
00084   fTreeOut->Branch("enu", &fENu, "enu/F");
00085   fTreeOut->Branch("emu", &fEMu, "mnu/F");
00086 
00087   fTreeOut->Branch("nret", &fNRetained, "nret/I");
00088   fTreeOut->Branch("nretmu", &fNRetainedMuon, "nretmu/I");
00089   fTreeOut->Branch("nretshw", &fNRetainedShw, "nretshw/I");
00090   fTreeOut->Branch("nretboth", &fNRetainedBoth, "nretboth/I");
00091 
00092   fTreeOut->Branch("peret", &fPERetained, "peret/F");
00093   fTreeOut->Branch("peretmu", &fPERetainedMuon, "peretmu/F");
00094   fTreeOut->Branch("peretshw", &fPERetainedShw, "peretshw/F");
00095   fTreeOut->Branch("peretboth", &fPERetainedBoth, "peretboth/F");
00096 
00097   fTreeOut->Branch("nrej", &fNRejected, "nrej/I");
00098   fTreeOut->Branch("nrejmu", &fNRejectedMuon, "nrejmu/I");
00099   fTreeOut->Branch("nrejshw", &fNRejectedShw, "nrejshw/I");
00100   fTreeOut->Branch("nrejboth", &fNRejectedBoth, "nrejboth/I");
00101 
00102   fTreeOut->Branch("nmuondig", &fNMuonDig, "nmuondig/I");
00103   fTreeOut->Branch("nmuondigret", &fNMuonDigRetained, "nmuondigret/I");
00104 
00105   fTreeOut->Branch("nshwdig", &fNShwDig, "nshwdig/I");
00106   fTreeOut->Branch("nshwdigret", &fNShwDigRetained, "nshwdigret/I");
00107 
00108   fTreeOut->Branch("nshwdigatvtx", &fNShwDigAtVtx, "nshwdigatvtx/I");
00109   fTreeOut->Branch("nshwdigretatvtx", &fNShwDigRetainedAtVtx, "nshwdigretatvtx/I");
00110 
00111   fTreeOut->Branch("nshwpe", &fNShwPE, "nshwpe/F");
00112   fTreeOut->Branch("nshwperet", &fNShwPERetained, "nshwpepe/F");
00113 
00114   fTreeOut->Branch("nshwpeatvtx", &fNShwPEAtVtx, "nshwpeatvtx/F");
00115   fTreeOut->Branch("nshwperetatvtx", &fNShwPERetainedAtVtx, "nshwpepeatvtx/F");
00116 
00117   fTreeOut->Branch("rshwn", &fRShwN, "rshwn/I");
00118   fTreeOut->Branch("rshwe", fRShwE, "rshwe[rshwn]/F");
00119   fTreeOut->Branch("rshwdist", fRShwDist, "rshwdist[rshwn]/I");
00120   fTreeOut->Branch("rshwmuon", fRShwMuon, "rshwmuon[rshwn]/I");
00121   fTreeOut->Branch("rshwtrk", fRShwTrk, "rshwtrk[rshwn]/I");
00122   fTreeOut->Branch("rshwxtalk", fRShwXTalk, "rshwxtalk[rshwn]/I");
00123 
00124   fTreeOut->Branch("rshwhittype", fRShwHitType, "rshwhittype[rshwn]/I");
00125   fTreeOut->Branch("rshwhitstatus", fRShwHitStatus, "rshwhitstatus[rshwn]/I");
00126   fTreeOut->Branch("rshwhitmother", fRShwHitMother, "rshwhitmoher[rshwn]/I");
00127 
00128   fTreeOut->Branch("fTrkMIPCalibOkay", &fTrkMIPCalibOkay, "fTrkMIPCalibOkay/I");
00129 
00130   hMIPTrueMuon  = new TH1F("hMIPTrueMuon", "", 100, 0, 10);
00131   hMIPMissTrack = new TH1F("hMIPMissTrack", "", 100, 0, 10);
00132   hMIPMixed     = new TH1F("hMIPMixed", "", 100, 0, 10);
00133   hMIPTrueMuon2  = new TH1F("hMIPTrueMuon2", "", 100, 0, 10);
00134   hMIPMissTrack2 = new TH1F("hMIPMissTrack2", "", 100, 0, 10);
00135   hMIPMixed2     = new TH1F("hMIPMixed2", "", 100, 0, 10);
00136  
00137   cdir->cd();
00138 }

void MCAnalysis::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 581 of file MCAnalysis.cxx.

References fMakePictures, fNameListIn, Registry::Get(), and Registry::GetCharString().

00582 {
00586   Int_t tmpi = 0;
00587   fNameListIn = r.GetCharString("NameListIn");
00588   if(r.Get("Draw", tmpi) ){
00589     fMakePictures = tmpi;
00590   }
00591 }

const Registry & MCAnalysis::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 558 of file MCAnalysis.cxx.

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

00559 {
00563   static Registry r; // Default configuration for module
00564 
00565   // Set name of config
00566   std::string name = this->GetName();
00567   name += ".config.default";
00568   r.SetName(name.c_str());
00569 
00570   // Set values in configuration
00571   r.UnLockValues();
00572   r.Set("NameListIn","mergedigitlist");
00573   r.Set("Draw",0);
00574   r.LockValues();
00575 
00576   return r;
00577 }

void MCAnalysis::EndJob (  )  [virtual]

Implement for notification of end of job

Reimplemented from JobCModule.

Definition at line 142 of file MCAnalysis.cxx.

References fFileOut, fTreeOut, hMIPMissTrack, hMIPMissTrack2, hMIPMixed, hMIPMixed2, hMIPTrueMuon, and hMIPTrueMuon2.

00143 {
00144   TDirectory* cdir = gDirectory;
00145   fFileOut->cd();
00146   fTreeOut->SetDirectory(fFileOut);
00147   fTreeOut->Write();
00148   hMIPTrueMuon ->SetDirectory(fFileOut);
00149   hMIPMissTrack->SetDirectory(fFileOut);
00150   hMIPMixed    ->SetDirectory(fFileOut);
00151   hMIPTrueMuon2->SetDirectory(fFileOut); 
00152   hMIPMissTrack2->SetDirectory(fFileOut);
00153   hMIPMixed2    ->SetDirectory(fFileOut);
00154 
00155   hMIPTrueMuon ->Write();
00156   hMIPMissTrack->Write();
00157   hMIPMixed    ->Write();
00158   hMIPTrueMuon2 ->Write();
00159   hMIPMissTrack2->Write();
00160   hMIPMixed2    ->Write();
00161   
00162   cdir->cd();
00163   
00164 }

void MCAnalysis::FillPlnInfo ( const CandDeMuxDigitListHandle digitlist  )  [private]

Definition at line 636 of file MCAnalysis.cxx.

References digit(), fPlnQ, CandHandle::GetDaughterIterator(), CalDigitType::kPE, and CandDeMuxDigit::kXTalk.

Referenced by Ana().

00637 {
00638   for(int i=0; i<500; i++) fPlnQ[i] = 0;
00639   TIter digititer(digitlist->GetDaughterIterator());
00640   while (const CandDeMuxDigitHandle *digit = dynamic_cast<const CandDeMuxDigitHandle*>(digititer())) {
00641       if(!(digit->GetDeMuxDigitFlagWord()&CandDeMuxDigit::kXTalk)) {
00642         const int plane = digit->GetPlexSEIdAltL().GetPlane();
00643         if(plane>0 && plane<500) fPlnQ[plane]+=digit->GetCharge(CalDigitType::kPE);
00644         
00645       } 
00646   }
00647 }

void MCAnalysis::FillTrkInfo ( const CandTrackHandle track  )  [private]

Definition at line 594 of file MCAnalysis.cxx.

References fTrkMIPDcosz, fTrkQ, fTrkStrips, fUZTrkQ, fVZTrkQ, CandHandle::GetDaughterIterator(), CandRecoHandle::GetEndDirCosZ(), CandRecoHandle::GetEndPlane(), CandRecoHandle::GetPlaneCharge(), CandRecoHandle::GetStripCharge(), CandTrackHandle::GetU(), CandTrackHandle::GetV(), CandRecoHandle::GetVtxDirCosZ(), CandRecoHandle::GetVtxPlane(), CandTrackHandle::GetZ(), CandTrackHandle::IsTPosValid(), CalStripType::kMIP, CalDigitType::kPE, and PlaneView::kU.

Referenced by Ana().

00594                                                         {
00595   for(int i=0; i<500; i++){
00596     fTrkQ[i] = 0;
00597     fTrkMIPDcosz[i] = 0;    
00598   }
00599   fTrkStrips.clear();
00600 
00601   const double nplanes =  track->GetEndPlane() - track->GetVtxPlane();
00602   for(int ipln = track->GetVtxPlane(); ipln<=track->GetEndPlane(); ++ipln){
00603     double dcosz =track->GetVtxDirCosZ() + ((ipln - track->GetVtxPlane())/nplanes)*(track->GetEndDirCosZ() - track->GetVtxDirCosZ() );
00604     int ip1 =-1;
00605     int ip2 =-1;
00606     if(track->IsTPosValid(ipln-1) && track->IsTPosValid(ipln+1)){
00607       ip2 = ipln+1; ip1 = ipln-1;
00608     }else if(track->IsTPosValid(ipln) && track->IsTPosValid(ipln+1)){
00609       ip2 = ipln+1; ip1 = ipln;
00610     }else if(track->IsTPosValid(ipln) && track->IsTPosValid(ipln-1)){
00611       ip2 = ipln; ip1 = ipln-1;
00612     }
00613     if(ip1!=-1 && ip2!=-1){
00614       const double l = sqrt((track->GetU(ip2) -  track->GetU(ip1))*(track->GetU(ip2) -  track->GetU(ip1)) 
00615                             + (track->GetV(ip2) -  track->GetV(ip1))*(track->GetV(ip2) -  track->GetV(ip1)) 
00616                             + (track->GetZ(ip2) -  track->GetZ(ip1))*(track->GetZ(ip2) -  track->GetZ(ip1)) ); 
00617       dcosz = fabs(track->GetZ(ip2) -  track->GetZ(ip1) ) / l;
00618     }
00619     //    cout << " plane: "<< ipln << "   " <<track->GetPlaneCharge(ipln, CalStripType::kMIP) <<endl; 
00620     fTrkMIPDcosz[ipln] = track->GetPlaneCharge(ipln, CalStripType::kMIP)*dcosz;
00621   }
00622 
00623 
00624   TIter trkstpIter(track->GetDaughterIterator());
00625   while(const CandStripHandle* strip = dynamic_cast<const CandStripHandle*>(trkstpIter())){
00626     fTrkStrips.insert(strip->GetPlane()*192+strip->GetStrip());
00627     fTrkQ[strip->GetPlane()]+=strip->GetCharge(CalDigitType::kPE);
00628     if(strip->GetPlaneView()==PlaneView::kU){
00629       fUZTrkQ->Fill(strip->GetPlane(), track->GetStripCharge(strip, CalStripType::kMIP) );
00630     }else{
00631       fVZTrkQ->Fill(strip->GetPlane(), track->GetStripCharge(strip, CalStripType::kMIP) );
00632     }
00633   }
00634 }

int MCAnalysis::GetMaxTrkPlane (  )  [private]

Definition at line 649 of file MCAnalysis.cxx.

References fPlnQ, and fTrkQ.

Referenced by Ana().

00649                               {
00650   int ntracklikeplanes=0;
00651   int maxtracklikeplane = 0;  
00652   for(int ipln = 0; ipln<500  && ntracklikeplanes<6; ipln++){
00653     if(fTrkQ[ipln]>0){
00654       if(fTrkQ[ipln]/fPlnQ[ipln]>.8 && fTrkQ[ipln]>3.0){
00655         ntracklikeplanes++;
00656         maxtracklikeplane=ipln;
00657       }
00658     }
00659   }
00660   if(ntracklikeplanes!=6) maxtracklikeplane = 500;
00661   return maxtracklikeplane;
00662 }


Member Data Documentation

Int_t MCAnalysis::fAction [private]

Definition at line 58 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

TCanvas* MCAnalysis::fCanvas [private]

Definition at line 111 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

TCanvas* MCAnalysis::fCanvas2 [private]

Definition at line 116 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fEMu [private]

Definition at line 60 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fENu [private]

Definition at line 59 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

TFile* MCAnalysis::fFileOut [private]

Definition at line 52 of file MCAnalysis.h.

Referenced by BeginJob(), and EndJob().

Int_t MCAnalysis::fMakePictures [private]

Definition at line 48 of file MCAnalysis.h.

Referenced by Ana(), and Config().

std::string MCAnalysis::fNameListIn [private]

Definition at line 47 of file MCAnalysis.h.

Referenced by Config().

Int_t MCAnalysis::fNMuonDig [private]

Definition at line 62 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Definition at line 63 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNRejected [private]

Definition at line 85 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNRejectedBoth [private]

Definition at line 88 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNRejectedMuon [private]

Definition at line 86 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNRejectedShw [private]

Definition at line 87 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNRetained [private]

Definition at line 75 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNRetainedBoth [private]

Definition at line 78 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNRetainedMuon [private]

Definition at line 76 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNRetainedShw [private]

Definition at line 77 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNShwDig [private]

Definition at line 65 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNShwDigAtVtx [private]

Definition at line 67 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Definition at line 66 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Definition at line 68 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fNShwPE [private]

Definition at line 70 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fNShwPEAtVtx [private]

Definition at line 72 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fNShwPERetained [private]

Definition at line 71 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Definition at line 73 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fNu [private]

Definition at line 57 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fPERetained [private]

Definition at line 80 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fPERetainedBoth [private]

Definition at line 83 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fPERetainedMuon [private]

Definition at line 81 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fPERetainedShw [private]

Definition at line 82 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fPlnQ[500] [private]

Definition at line 126 of file MCAnalysis.h.

Referenced by FillPlnInfo(), and GetMaxTrkPlane().

Int_t MCAnalysis::fRShwDist[500] [private]

Definition at line 92 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fRShwE[500] [private]

Definition at line 91 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fRShwHitMother[500] [private]

Definition at line 99 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fRShwHitStatus[500] [private]

Definition at line 98 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fRShwHitType[500] [private]

Definition at line 97 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fRShwMuon[500] [private]

Definition at line 93 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fRShwN [private]

Definition at line 90 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fRShwTrk[500] [private]

Definition at line 94 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fRShwXTalk[500] [private]

Definition at line 95 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fRun [private]

Definition at line 55 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Int_t MCAnalysis::fSnarl [private]

Definition at line 56 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

TTree* MCAnalysis::fTreeOut [private]

Definition at line 53 of file MCAnalysis.h.

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

Definition at line 100 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

Float_t MCAnalysis::fTrkMIPDcosz[500] [private]

Definition at line 125 of file MCAnalysis.h.

Referenced by Ana(), and FillTrkInfo().

Float_t MCAnalysis::fTrkQ[500] [private]

Definition at line 124 of file MCAnalysis.h.

Referenced by FillTrkInfo(), and GetMaxTrkPlane().

set<int> MCAnalysis::fTrkStrips [private]

Definition at line 130 of file MCAnalysis.h.

Referenced by Ana(), and FillTrkInfo().

vector<TMarker> MCAnalysis::fUHits [private]

Definition at line 114 of file MCAnalysis.h.

Referenced by Ana().

TH2F* MCAnalysis::fUZAxis [private]

Definition at line 112 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

TH1F* MCAnalysis::fUZTrkQ [private]

Definition at line 117 of file MCAnalysis.h.

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

vector<TMarker> MCAnalysis::fVHits [private]

Definition at line 115 of file MCAnalysis.h.

Referenced by Ana().

TH2F* MCAnalysis::fVZAxis [private]

Definition at line 113 of file MCAnalysis.h.

Referenced by Ana(), and BeginJob().

TH1F* MCAnalysis::fVZTrkQ [private]

Definition at line 118 of file MCAnalysis.h.

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

TH1F* MCAnalysis::hMIPMissTrack [private]

Definition at line 103 of file MCAnalysis.h.

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

TH1F* MCAnalysis::hMIPMissTrack2 [private]

Definition at line 106 of file MCAnalysis.h.

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

TH1F* MCAnalysis::hMIPMixed [private]

Definition at line 104 of file MCAnalysis.h.

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

TH1F* MCAnalysis::hMIPMixed2 [private]

Definition at line 107 of file MCAnalysis.h.

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

TH1F* MCAnalysis::hMIPTrueMuon [private]

Definition at line 102 of file MCAnalysis.h.

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

TH1F* MCAnalysis::hMIPTrueMuon2 [private]

Definition at line 105 of file MCAnalysis.h.

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


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

Generated on 2 Nov 2017 for loon by  doxygen 1.6.1