Inheritance diagram for icarus::crt::CRTTrackProducer:

Public Member Functions
	CRTTrackProducer (fhicl::ParameterSet const &p)

	CRTTrackProducer (CRTTrackProducer const &)=delete

	CRTTrackProducer (CRTTrackProducer &&)=delete

CRTTrackProducer &	operator= (CRTTrackProducer const &)=delete

CRTTrackProducer &	operator= (CRTTrackProducer &&)=delete

void	produce (art::Event &e) override

void	beginJob () override

void	endJob () override

Private Attributes
string	fDataLabelHits

int	fTrackMethodType

int	fStoreTrack

bool	fUseTopPlane

CRTTrackRecoAlg	trackAlg

Detailed Description

Definition at line 50 of file icaruscode/icaruscode/CRT/CRTTrackProducer_module.cc.

Constructor & Destructor Documentation

icarus::crt::CRTTrackProducer::CRTTrackProducer ( fhicl::ParameterSet const & p )

explicit

Definition at line 115 of file icaruscode/icaruscode/CRT/CRTTrackProducer_module.cc.

     : EDProducer{p}, trackAlg(p.get<fhicl::ParameterSet>("TrackAlg"))
   {  
       // Initialize member data here.
       fDataLabelHits      = p.get<string>("DataLabelHits");      // CRTHit producer module name
       // fDataLabelTZeros    = p.get<string>("DataLabelTZeros");    // CRTTzero producer module name
       fStoreTrack         = p.get<int>   ("StoreTrack");         // method 1 = all, method 2 = ave, method 3 = pure, method 4 = top plane
       fTrackMethodType    = p.get<int>   ("TrackMethodType");    // Print stuff
       fUseTopPlane        = p.get<bool>  ("UseTopPlane");        // Use hits from the top plane (SBND specific)
       
       // Call appropriate produces<>() functions here.
       if(fStoreTrack == 1){ 
           produces< vector<sbn::crt::CRTTrack>   >();
           produces< art::Assns<sbn::crt::CRTTrack , sbn::crt::CRTHit> >();
       } 
     
   } // CRTTrackProducer()

icarus::crt::CRTTrackProducer::CRTTrackProducer ( CRTTrackProducer const & )

delete

icarus::crt::CRTTrackProducer::CRTTrackProducer ( CRTTrackProducer && )

delete

Member Function Documentation

void icarus::crt::CRTTrackProducer::beginJob ( )

override

Definition at line 385 of file icaruscode/icaruscode/CRT/CRTTrackProducer_module.cc.

386 {

387 } // CRTTrackProducer::beginJob()

void icarus::crt::CRTTrackProducer::endJob ( )

override

Definition at line 389 of file icaruscode/icaruscode/CRT/CRTTrackProducer_module.cc.

390 {

391 } // CRTTrackProducer::endJob()

CRTTrackProducer& icarus::crt::CRTTrackProducer::operator= ( CRTTrackProducer const & )

delete

CRTTrackProducer& icarus::crt::CRTTrackProducer::operator= ( CRTTrackProducer && )

delete

void icarus::crt::CRTTrackProducer::produce ( art::Event & e )

override

Definition at line 133 of file icaruscode/icaruscode/CRT/CRTTrackProducer_module.cc.

   {
       // For validation
       int nTrack = 0;
       int nCompTrack = 0;
       int nIncTrack = 0;
     
       // CRTTrack collection on this event                                                                         
       std::unique_ptr<vector<sbn::crt::CRTTrack> > CRTTrackCol(new vector<sbn::crt::CRTTrack>);
       std::unique_ptr< art::Assns<sbn::crt::CRTTrack, sbn::crt::CRTHit> > Trackassn( new art::Assns<sbn::crt::CRTTrack, sbn::crt::CRTHit>);
       art::PtrMaker<sbn::crt::CRTTrack> makeTrackPtr(evt);
     
       // Implementation of required member function here.
       art::Handle< vector<sbn::crt::CRTHit> > rawHandle;
       evt.getByLabel(fDataLabelHits, rawHandle); //what is the product instance name? no BernZMQ
     
       // Check to make sure the data we asked for is valid                                                                         
       if(!rawHandle.isValid()){
           mf::LogWarning("CRTTrackProducer")
             <<"No CRTHits from producer module "<<fDataLabelHits;
           return;
       }
     
       // Track method 4 = SBND method with top plane (doesn't use CRTTzero)
       if(fTrackMethodType == 4){
     
           //Get the CRT hits from the event
           vector<art::Ptr<sbn::crt::CRTHit> > hitlist;
     
           if (evt.getByLabel(fDataLabelHits, rawHandle))
             art::fill_ptr_vector(hitlist, rawHandle);
     
           map<art::Ptr<sbn::crt::CRTHit>, int> hitIds;
     
           for(size_t i = 0; i<hitlist.size(); i++){
               hitIds[hitlist[i]] = i;
           }
     
           vector<vector<art::Ptr<sbn::crt::CRTHit>>> CRTTzeroVect = trackAlg.CreateCRTTzeros(hitlist);
     
           // Loop over tzeros
           for(size_t i = 0; i<CRTTzeroVect.size(); i++){
     
               //loop over hits for this tzero, sort by tagger
               //map<string, vector<art::Ptr<sbn::crt::CRTHit>>> hits;
             map<int, vector<art::Ptr<sbn::crt::CRTHit>>> hits;
     
               for (size_t ah = 0; ah< CRTTzeroVect[i].size(); ++ah){        
     
                 //string ip = CRTTzeroVect[i][ah]->tagger;       
                 int ip = CRTTzeroVect[i][ah]->plane;       
                 hits[ip].push_back(CRTTzeroVect[i][ah]);
               } // loop over hits
               
               //loop over planes and calculate average hits
               vector<pair<sbn::crt::CRTHit, vector<int>>> allHits;
     
               for (auto &keyVal : hits){
     
                 // string ip = keyVal.first;
                 int ip = keyVal.first;
                 vector<pair<sbn::crt::CRTHit, vector<int>>> ahits = trackAlg.AverageHits(hits[ip], hitIds);
 
                 allHits.insert(allHits.end(), ahits.begin(), ahits.end());
     
                 /*
                   if(fUseTopPlane && ip == "volTaggerTopHigh_0"){ 
                       allHits.insert(allHits.end(), ahits.begin(), ahits.end());
                   }
     
                   else if(ip != "volTaggerTopHigh_0"){ 
                       allHits.insert(allHits.end(), ahits.begin(), ahits.end());
                   }
                 */
               }
     
               //Create tracks with hits at the same tzero
               vector<pair<sbn::crt::CRTTrack, vector<int>>> trackCandidates = trackAlg.CreateTracks(allHits);
               nTrack += trackCandidates.size();
     
               for(size_t j = 0; j < trackCandidates.size(); j++){
     
                   CRTTrackCol->emplace_back(trackCandidates[j].first);
                   art::Ptr<sbn::crt::CRTTrack> trackPtr = makeTrackPtr(CRTTrackCol->size()-1);
     
                   for (size_t ah = 0; ah< CRTTzeroVect[i].size(); ++ah){        
                       Trackassn->addSingle(trackPtr, CRTTzeroVect[i][ah]);
                   }
     
                   if(trackCandidates[j].first.complete) 
                       nCompTrack++;
                   else 
                       nIncTrack++;
               }
           }
     
       }//end if sbnd option
 
        /*   
       //Older track reconstruction methods from MicroBooNE
       else{
           //Get list of tzeros             
           //  std::vector<crt::CRTHit> const& CRTHitCollection(*rawHandle);
           art::Handle< vector<sbn::crt::CRTTzero> > rawHandletzero;
           evt.getByLabel(fDataLabelTZeros, rawHandletzero); //what is the product instance name? no BernZMQ
           
           //check to make sure the data we asked for is valid                                                                            
           if(!rawHandletzero.isValid()){
             mf::LogWarning("CRTTrackProducer")
               <<"No CRTTzeros from producer module "<<fDataLabelTZeros;
             return;
           }
           
           vector<art::Ptr<sbn::crt::CRTTzero> > tzerolist;
           if (evt.getByLabel(fDataLabelTZeros,rawHandletzero))
             art::fill_ptr_vector(tzerolist, rawHandletzero);
      
           art::FindManyP<sbn::crt::CRTHit> fmht(rawHandletzero, evt, fDataLabelTZeros);
      
           //loop over tzeros
           for(size_t tzIter = 0; tzIter < tzerolist.size(); ++tzIter){   
             
               //count planes with hits for this tzero
               int np =0 ; //int ipflag[7] = {}; // CHANGED FROM 4 TO 7
               int tothits =0;
               for (int ip=0;ip<7;++ip) { // CHANGED FROM 4 TO 7
                   if (tzerolist[tzIter]->nhits[ip]>0){ 
                       np++; 
                       //ipflag[ip]=1; 
                       tothits+=tzerolist[tzIter]->nhits[ip];
                   }  
               }
      
               if (np<2) continue;
               vector<art::Ptr<sbn::crt::CRTHit> > hitlist=fmht.at(tzIter);
               //for(size_t hit_i = 0; hit_i < hitlist.size(); hit
               if (fTrackMethodType==1) {
                   double time_s_A = hitlist[0]->ts0_s;
      
                   // find pairs of hits in different planes
                   for (size_t ah = 0; ah< hitlist.size()-1; ++ah){        
                       sbn::crt::CRTHit temphit=*hitlist[ah];
                       CRTavehit Ahit = copyme(temphit);
                       int planeA = hitlist[ah]->plane;
      
                       for (size_t bh = ah+1; bh< hitlist.size(); ++bh){        
                           int planeB = hitlist[bh]->plane;
      
                           if (planeB!=planeA && !((planeA==3&&planeB==4)||(planeA==4&&planeB==3))) {  // make a track               
                               temphit=*hitlist[bh];
                               CRTavehit Bhit = copyme(temphit);
                               sbn::crt::CRTTrack CRTcanTrack=shcut(Ahit,Bhit,time_s_A,0);
                               CRTTrackCol->emplace_back(CRTcanTrack);
                           }
      
                       }
      
                   }
      
               }
      
               else if ((fTrackMethodType==2) || (fTrackMethodType==3 && np==2 && tothits==2)) {        
      
                   //loop over hits and get average x,y,z,pe for each plane CHANGED FROM 4 TO 7
                   vector<float> thittime0[7];
                   vector<float> thittime1[7];
                   vector<float> tx[7];
                   vector<float> ty[7];
                   vector<float> tz[7];
                   vector<float> pe[7];
                   
                   double time_s_A = hitlist[0]->ts0_s;
                   //      double time_s_err = hitlist[0]->ts0_s_err;
                   double time_s_err = 0.;
                   double time1_ns_A = hitlist[0]->ts1_ns;
                   double time0_ns_A = hitlist[0]->ts0_ns;
                
                   //loop over hits for this tzero, sort by plane
                   for (size_t ah = 0; ah< hitlist.size(); ++ah){        
                       int ip = hitlist[ah]->plane;       
                       thittime0[ip].push_back(hitlist[ah]->ts0_ns-time0_ns_A);
                       thittime1[ip].push_back(hitlist[ah]->ts1_ns-time1_ns_A);
                       tx[ip].push_back(hitlist[ah]->x_pos);
                       ty[ip].push_back(hitlist[ah]->y_pos);
                       tz[ip].push_back(hitlist[ah]->z_pos);
                       pe[ip].push_back(hitlist[ah]->peshit);        
                   } // loop over hits
      
                   CRTavehit aveHits[7];
                   //loop over planes and calculate average hits
                   for (int ip = 0; ip < 7; ip++){
                       if (tx[ip].size()>0){
                           uint32_t at0; int32_t at1; uint16_t rt0,rt1;
                           float totpe=0.0;
                           float avet1=0.0; float rmst1 =0.0; 
                           float avet0=0.0; float rmst0 =0.0; 
                           float avex=0.0; float rmsx =0.0; 
                           float avey=0.0; float rmsy =0.0; 
                           float avez=0.0; float rmsz =0.0;
                           vmanip(thittime0[ip],&avet0,&rmst0);
                           vmanip(thittime1[ip],&avet1,&rmst1);
                           at0 = (uint32_t)(avet0+time0_ns_A); rt0 = (uint16_t)rmst0;   
                           at1 = (int32_t)(avet1+time1_ns_A); rt1 = (uint16_t)rmst1;
                           vmanip(tx[ip],&avex,&rmsx);
                           vmanip(ty[ip],&avey,&rmsy);
                           vmanip(tz[ip],&avez,&rmsz);
                           totpe=std::accumulate(pe[ip].begin(), pe[ip].end(), 0.0);
                           CRTavehit aveHit = fillme(at0,rt0,at1,rt1,avex,rmsx,avey,rmsy,avez,rmsz,totpe,ip,"");
                           aveHits[ip] = aveHit;
                       }
                       else {
                           CRTavehit aveHit = fillme(0,0,0,0,-99999,-99999,-99999,-99999,-99999,-99999,-99999,ip,"");
                           aveHits[ip] = aveHit;
                       }
                   }  
      
                   // find pairs of hits in different planes
                   for (size_t ah = 0; ah< 6; ++ah){        
                       CRTavehit Ahit = aveHits[ah];
                       if( Ahit.x_pos==-99999 ) continue;
      
                       for (size_t bh = ah+1; bh< 7; ++bh){        
                           if ( aveHits[bh].x_pos==-99999 ) continue;
      
                           if (ah!=bh && !(ah==3&&bh==4)) {  // make a track               
                               CRTavehit Bhit = aveHits[bh];
                               sbn::crt::CRTTrack CRTcanTrack=shcut(Ahit,Bhit,time_s_A,time_s_err);
                               CRTTrackCol->emplace_back(CRTcanTrack);
                               nTrack++;
                           }
      
                       }
      
                   }
      
               }
             
           }// loop over tzeros
       }//uBooNE method
       */
       //store track collection into event
       if(fStoreTrack == 1){
           evt.put(std::move(CRTTrackCol));
           evt.put(std::move(Trackassn));
       }
       mf::LogInfo("CRTTrackProducer")
         <<"Number of tracks            = "<<nTrack<<"\n"
         <<"Number of complete tracks   = "<<nCompTrack<<"\n"
         <<"Number of incomplete tracks = "<<nIncTrack;
       
   } // CRTTrackProducer::produce()