#include <ElecSimProc.h>
Inheritance diagram for ElecSimProc:
Based on EsFrontEndAlg.
This algorithms handles the care and feeding of ElecSim in a staged processing ("pull") context.
ElecSim requires SimHitCollections to be suitably chunked in time. Chunks are built by collecting hits ordered in time until a gap of at least "preTimeTolerance + postTimeTolerance" can be found in the hits from all detectors in a site. What chunks have been built are packaged in to a fake SimHitHeader and sent through ElecSim for to-pulse and to-crate processing. The results are packaged into a ElecHeader. Any unused SimHits are cached for subsequent use.
bv@bnl.gov Thu Oct 30 17:05:59 2008
Here hits are made into some sets where there is a gap between them. The hits from one SimHeader can split into several sets and the hits from different SimHeaders can overlap.
Two parameters are added. One is chunk start time. The other one is chunk stop time. They are used to trim SimHit pipeline and split the hits in one SimHeader. Like IBD event might be split into to two fake SimHitHeader.
Even after a time gap is found, this algorithm still need to keep pulling SimHeader from lower stage Detector. Like for IBD event, after a time gap is found there is a chance the hits from next SimHeader is going to overlay with the prompt positron signal.
The approach of ElecSimProc is different with DetSimProc Since a time gap arount 10us is needed to find between two sets of successive hits, then don't have to worry about the time crossing caused by electronic simulation.
Zhe Wang rewritten at Dec 21 2008
MuonProphet can produce some fast simulation result for muons. Empty ElecHeaders are generated with inputHeader pointing back to the GenHeader of these muons. These empty ElecHeaders are passed to trigger and readout simulation right after they are created. Those fast simulated muon are called noble particle with respect to Geant4.
Zhe Wang Mar 11 2010
Definition at line 74 of file ElecSimProc.h.
typedef HeaderStageData<DayaBay::HeaderObject> ElecSimProc::HeaderData [private] |
Definition at line 116 of file ElecSimProc.h.
typedef HeaderStageData<DayaBay::SimHeader> ElecSimProc::SimData [private] |
Definition at line 117 of file ElecSimProc.h.
typedef HeaderStageData<DayaBay::ElecHeader> ElecSimProc::ElecData [private] |
Definition at line 118 of file ElecSimProc.h.
typedef std::multimap<TimeStamp,DayaBay::SimHit*, std::less<TimeStamp> > ElecSimProc::HitPipeline [private] |
Definition at line 173 of file ElecSimProc.h.
ElecSimProc::ElecSimProc | ( | const std::string & | name, | |
ISvcLocator * | pSvcLocator | |||
) |
Definition at line 36 of file ElecSimProc.cc.
00037 : StageProcessor<DayaBay::ElecHeader>(name,pSvcLocator) 00038 , m_pmtTool(0) 00039 , m_rpcTool(0) 00040 , m_feeTool(0) 00041 , m_fecTool(0) 00042 { 00043 00044 m_pullMode=true; 00045 00046 // ElecSim 00047 declareProperty("Detectors",m_detectorNames,"List of active detectors"); 00048 declareProperty("PmtTool",m_pmtToolName="EsPmtEffectPulseTool", 00049 "Name of the PMT simulation tool"); 00050 declareProperty("RpcTool",m_rpcToolName="EsIdealPulseTool", 00051 "Name of the RPC simulation tool"); 00052 declareProperty("FeeTool",m_feeToolName="EsIdealFeeTool", 00053 "Name of the PMT Front-end electronics simulation tool"); 00054 declareProperty("FecTool",m_fecToolName="EsIdealFecTool", 00055 "Name of the RPC Front-end electronics simulation tool"); 00056 // Default configuration is DayaBay near site, ADs only. 00057 m_detectorNames.push_back("DayaBayAD1"); 00058 m_detectorNames.push_back("DayaBayAD2"); 00059 m_detectorNames.push_back("DayaBayIWS"); 00060 m_detectorNames.push_back("DayaBayOWS"); 00061 m_detectorNames.push_back("DayaBayRPC"); 00062 m_detectorNames.push_back("LingAoAD1"); 00063 m_detectorNames.push_back("LingAoAD2"); 00064 m_detectorNames.push_back("LingAoIWS"); 00065 m_detectorNames.push_back("LingAoOWS"); 00066 m_detectorNames.push_back("LingAoRPC"); 00067 m_detectorNames.push_back("FarAD1"); 00068 m_detectorNames.push_back("FarAD2"); 00069 m_detectorNames.push_back("FarAD3"); 00070 m_detectorNames.push_back("FarAD4"); 00071 m_detectorNames.push_back("FarIWS"); 00072 m_detectorNames.push_back("FarOWS"); 00073 m_detectorNames.push_back("FarRPC"); 00074 00075 // special for this package 00076 declareProperty("MinGapTime", m_minTimeGap = 00077 DayaBay::preTimeTolerance + DayaBay::postTimeTolerance, 00078 "Minimum width for a gap that defines a chunk of data."); 00079 declareProperty("LongestJump", m_longest = 10*365*24*60*60, 00080 "A protection for longest jump in time stamp. Default is set to 10 year."); 00081 00082 m_hitPipeline.clear(); 00083 m_giantHits.clear(); 00084 m_giantHitsInTime.clear(); 00085 m_giantHitSheader.clear(); 00086 00087 m_currentTime = TimeStamp::GetBOT(); 00088 m_realElecHeaderTime = TimeStamp::GetBOT(); 00089 m_emptyElecHeaderTime = TimeStamp::GetBOT(); 00090 }
ElecSimProc::~ElecSimProc | ( | ) | [virtual] |
StatusCode ElecSimProc::initialize | ( | ) | [virtual] |
Reimplemented from StageProcessor< DayaBay::ElecHeader >.
Definition at line 96 of file ElecSimProc.cc.
00097 { 00098 StatusCode sc = this->StageProcessor<DayaBay::ElecHeader>::initialize(); 00099 if (sc.isFailure()) return sc; 00100 00102 // Convert detector names to Detector IDs 00103 int nDetNames = m_detectorNames.size(); 00104 for(int detIdx=0; detIdx<nDetNames; detIdx++){ 00105 DayaBay::Detector det(m_detectorNames[detIdx]); 00106 if(det.site() == Site::kUnknown 00107 || det.detectorId() == DetectorId::kUnknown){ 00108 error() << "Invalid detector name: " << m_detectorNames[detIdx] << endreq; 00109 return StatusCode::FAILURE; 00110 } 00111 m_detectors.push_back(det); 00112 } 00113 00114 // Get PMT Pulse tool 00115 try { 00116 m_pmtTool = tool<IEsPulseTool>(m_pmtToolName); 00117 } 00118 catch(const GaudiException& exg) { 00119 fatal() << "Failed to get pmt tool: \"" << m_pmtToolName << "\"" << endreq; 00120 return StatusCode::FAILURE; 00121 } 00122 info () << "Added tool " << m_pmtToolName << endreq; 00123 00124 // Get RPC Pulse tool 00125 if(m_pmtToolName == m_rpcToolName){ 00126 m_rpcTool = m_pmtTool; 00127 }else { 00128 try { 00129 m_rpcTool = tool<IEsPulseTool>(m_rpcToolName); 00130 } 00131 catch(const GaudiException& exg) { 00132 fatal() << "Failed to get rpc tool: \"" << m_rpcToolName << "\"" << endreq; 00133 return StatusCode::FAILURE; 00134 } 00135 info () << "Added tool " << m_rpcToolName << endreq; 00136 } 00137 00138 // Get PMT Front-end Electronics tool 00139 try { 00140 m_feeTool = tool<IEsFrontEndTool>(m_feeToolName); 00141 } 00142 catch(const GaudiException& exg) { 00143 fatal() << "Failed to get fee tool: \"" << m_feeToolName << "\"" << endreq; 00144 return StatusCode::FAILURE; 00145 } 00146 info () << "Added tool " << m_feeToolName << endreq; 00147 00148 // Get RPC Front-end Electronics tool 00149 try { 00150 m_fecTool = tool<IEsFrontEndTool>(m_fecToolName); 00151 } 00152 catch(const GaudiException& exg) { 00153 fatal() << "Failed to get fec tool: \"" << m_fecToolName << "\"" << endreq; 00154 return StatusCode::FAILURE; 00155 } 00156 info () << "Added tool " << m_fecToolName << endreq; 00157 00158 // Set the 'partition' for this run in RunDataSvc 00159 IRunDataSvc* runDataSvc = svc<IRunDataSvc>("RunDataSvc",true); 00160 if(!runDataSvc){ 00161 error() << "Failed to get RunDataSvc" << endreq; 00162 return StatusCode::FAILURE; 00163 } 00164 int task = 1; // Tells RunDataSvc to use new simulation run data. 00165 Context emptyContext; // Need a dummy context for this service call 00166 ServiceMode svcMode(emptyContext, task); 00167 const DayaBay::RunData* runData = runDataSvc->runData(svcMode); 00168 if(!runData){ 00169 error() << "Failed to get RunData for simulation" << endreq; 00170 return StatusCode::FAILURE; 00171 } 00172 DayaBay::RunData modifiedRunData(*runData); 00173 modifiedRunData.setDetectors(m_detectors); 00174 sc = runDataSvc->setRunData(modifiedRunData); 00175 if( !sc.isSuccess() ){ 00176 error() << "Failed to set detectors in run" << endreq; 00177 } 00178 00179 info()<<"minTimeGap "<<m_minTimeGap << " ns" <<endreq; 00180 00181 return sc; 00182 }
StatusCode ElecSimProc::execute | ( | ) | [virtual] |
Definition at line 184 of file ElecSimProc.cc.
00185 { 00186 debug()<<"Executing ..."<<endreq; 00187 // Only do anything if the current stage time has advanced beyond 00188 // the last time we ran. 00189 static bool first_time = true; 00190 if (!first_time && m_currentTime > thisStage()->currentTime()) { 00191 return StatusCode::SUCCESS; 00192 } 00193 first_time = false; 00194 00195 StatusCode sc; 00196 00197 sc = this->findGap(); 00198 if (sc.isFailure()) return sc; 00199 00200 sc = this->fillHitHeader(); 00201 if (sc.isFailure()) return sc; 00202 00203 sc = this->findGiantHitInTime(); 00204 if (sc.isFailure()) return sc; 00205 00206 // Number of normal hits is above zero 00207 if( m_count>0 ) { 00208 sc = this->preExecute(); 00209 if (sc.isFailure()) return sc; 00210 sc = this->ES_execute(); 00211 if (sc.isFailure()) return sc; 00212 sc = this->postExecute(); 00213 if (sc.isFailure()) return sc; 00214 } 00215 00216 // Number of giant hits is above zero 00217 if( m_giantHitsInTime.size()>0 ) { 00218 sc = this->fastElecSim(); 00219 if (sc.isFailure()) return sc; 00220 } 00221 00222 // update m_currentTime 00223 if( m_realElecHeaderTime > m_emptyElecHeaderTime ) { 00224 m_currentTime = m_realElecHeaderTime; 00225 } else { 00226 m_currentTime = m_emptyElecHeaderTime; 00227 } 00228 00229 // blowChunks() must be the last one, because it will hold the last refcount of some 00230 // input headers. Before do this, it must be added to inputheaders in postExecute(). 00231 sc = this->blowChunks(); 00232 if (sc.isFailure()) return sc; 00233 00234 return StatusCode::SUCCESS; 00235 }
StatusCode ElecSimProc::getSimData | ( | DayaBay::SimHeader *& | output | ) | [private] |
Definition at line 238 of file ElecSimProc.cc.
00239 { 00240 // Currently we get it from the lower stage, but this function 00241 // will eventually handle getting it from AES when the SimHeaders 00242 // are read in from file. 00243 00244 IStageData* stageData = 0; 00245 StatusCode sc = lowerStage()->nextElement(stageData); 00246 if (sc.isFailure()) { 00247 error() << "Failed to get data from detector simulation stage" 00248 << endreq; 00249 return sc; 00250 } 00251 00252 DayaBay::SimHeader* shead = 0; 00253 00255 SimData* input=0; 00256 input = dynamic_cast<SimData*>(stageData); 00257 00258 if( input!=0 ) { 00259 shead = &input->header(); 00260 } else { 00262 HeaderData* pHeaderData = dynamic_cast<HeaderData*>(stageData); 00263 DayaBay::HeaderObject* pHeaderObj = &pHeaderData->header(); 00264 shead = dynamic_cast<DayaBay::SimHeader*>(pHeaderObj); 00265 } 00266 00267 if( shead==0 ) { 00268 error() << "Failed to get SimHeader from lower stage" <<endreq; 00269 } 00270 00271 debug() << "stage data's time "<< stageData->time()<<endreq; 00272 00273 // we must temporarily hold on to this DataObject because this 00274 // input may have been created during previous exection cycle. 00275 shead->addRef(); // And because deleting the SimData 00276 delete stageData; // drops the only other reference. It 00277 stageData=0; // will be released once it is no 00278 // longer providing any hits to the 00279 // pipeline 00280 00281 debug() << "Current time of lower stage "<<lowerStage()->currentTime()<<endreq; 00282 debug() << "Got SimHeader from lower stage at " << shead->timeStamp() << endreq; 00283 00284 output = shead; 00285 return StatusCode::SUCCESS; 00286 }
StatusCode ElecSimProc::fillPipeLine | ( | DayaBay::SimHeader * | shead | ) | [private] |
Definition at line 290 of file ElecSimProc.cc.
00291 { 00292 // pipeline 00293 TimeStamp realLatest = TimeStamp::GetBOT(); 00294 00295 // Get the SimHit Header 00296 const DayaBay::SimHitHeader* shitHeader = shead->hits(); 00297 if (!shitHeader) { 00298 warning() << "Got NULL SimHitHeader with SimHeader @ " 00299 << shead->timeStamp() 00300 << " skipping..." 00301 << endreq; 00302 return StatusCode::SUCCESS; 00303 } 00304 const DayaBay::SimHitHeader::hc_map& hcMap = shitHeader->hitCollection(); 00305 00306 const TimeStamp& absTime = shead->timeStamp(); 00307 debug()<<"absTime = "<<absTime<<endreq; 00309 DayaBay::SimHitHeader::hc_map::const_iterator hcIter; 00310 for (hcIter=hcMap.begin(); hcIter != hcMap.end(); ++hcIter) { 00311 00312 DayaBay::Detector det(hcIter->first); 00313 DayaBay::SimHitCollection* hits = hcIter->second; 00314 if(!hits) return StatusCode::FAILURE; 00315 00316 const DayaBay::SimHitCollection::hit_container& hc = hits->collection(); 00317 for (size_t ind=0; ind < hc.size(); ++ind) { 00318 DayaBay::SimHit* hit = hc[ind]; 00319 00320 DayaBay::AdPmtSensor pmtid(hit->sensDetId()); 00321 if (pmtid.bogus()) { 00322 warning() << "SimHeaderCnvHits: warning got bogus PMTid: " << pmtid << endreq; 00323 } 00324 00325 TimeStamp ts = absTime; 00326 if( hit->hitTime()/CLHEP::second < m_longest ) { 00327 ts.Add(hit->hitTime()/CLHEP::second); 00328 m_hitPipeline.insert(pair<TimeStamp,DayaBay::SimHit*>(ts,hit)); 00329 //debug()<<hit->hitTime()/CLHEP::second<<" [sec]"<<endreq; 00330 // find the real latest time 00331 if( ts>realLatest ) realLatest = ts; 00332 } else { 00333 warning()<<"A hit in the far future "<<hit->hitTime()/CLHEP::second<<" [sec] is skipped."<<endreq; 00334 } 00335 } 00336 } 00337 00342 const DayaBay::IHeader* iheader = (shead->inputHeaders()) [0]; 00343 const DayaBay::GenHeader* genHeader = dynamic_cast<const DayaBay::GenHeader*>(iheader); 00344 const HepMC::GenEvent* genEvent = genHeader->event(); 00345 00346 if( genEvent ) { 00347 HepMC::GenEvent::particle_const_iterator pci, pci_end = genEvent->particles_end(); 00348 00349 for( pci = genEvent->particles_begin(); pci != pci_end; ++pci ) { 00350 if( (*pci)->pdg_id() == 13 || (*pci)->pdg_id() == -13 ) { // muon 00351 if( (*pci)->status() != MpMuonFate::kNeedSim ) { // geant4-noble particle 00352 debug()<<"Find one geant4-noble particle"<<endreq; 00354 DayaBay::SimHit* aGiantHit = new DayaBay::SimHit; 00357 double dt = 0 * CLHEP::nanosecond; 00358 aGiantHit->setHitTime(dt); 00360 aGiantHit->setWeight(0); 00361 aGiantHit->setSensDetId( DetectorId::kUnknown ); 00362 debug()<<"absTime "<<absTime<<endreq; 00363 TimeStamp ts = absTime; 00364 ts.Add(dt/CLHEP::second); 00365 debug()<<"ts "<<ts<<endreq; 00366 m_hitPipeline.insert(pair<TimeStamp,DayaBay::SimHit*>(ts,aGiantHit)); 00368 m_giantHits.insert(pair<TimeStamp,DayaBay::SimHit*>(ts,aGiantHit)); 00369 m_giantHitSheader[aGiantHit] = shead; 00370 debug()<<"aGiantHit is added "<<aGiantHit<<" with SimHeader "<<shead<<" at time "<<ts<<endreq; 00371 debug()<<"Now m_giantHitSheader has "<<m_giantHitSheader.size()<<" giant hits."<<endreq; 00372 // find the real latest time 00373 if( ts>realLatest ) realLatest = ts; 00374 break; // only process one fast simulated muon 00375 } 00376 } 00377 } 00378 } 00379 00380 debug()<<"Current m_giantHitSheader begin"<<endreq; 00381 map<const DayaBay::SimHit*, const DayaBay::SimHeader*>::iterator it,itend = m_giantHitSheader.end(); 00382 for( it=m_giantHitSheader.begin(); it!=itend; it++ ) { 00383 debug()<<"aGiantHit "<<it->first<<" with SimHeader "<<it->second<<endreq; 00384 } 00385 00386 // current SimHeaders 00387 m_currentHeaders.push_back(shead); 00388 debug()<<"m_currentHeaders has "<<m_currentHeaders.size()<<" SimHeaders."<<endreq; 00389 00390 // record the real latest time 00391 m_realLatest.insert( pair<DayaBay::SimHeader*, TimeStamp>(shead,realLatest) ); 00392 00393 return StatusCode::SUCCESS; 00394 }
StatusCode ElecSimProc::moar | ( | ) | [private] |
Definition at line 396 of file ElecSimProc.cc.
00397 { 00398 debug() << "Getting MOAR!" << endreq; 00399 00400 DayaBay::SimHeader* shead = 0; 00401 StatusCode sc = this->getSimData(shead); 00402 if (sc.isFailure()) return sc; 00403 00404 sc = this->fillPipeLine(shead); 00405 if (sc.isFailure()) return sc; 00406 00407 return StatusCode::SUCCESS; 00408 }
StatusCode ElecSimProc::findGap | ( | ) | [private] |
Definition at line 410 of file ElecSimProc.cc.
00411 { 00412 debug()<<"Entering findGap()"<<endreq; 00413 while ( m_hitPipeline.size() < 2 ) { // at least two hits are needed to find a gap 00414 StatusCode sc = this->moar(); // need MOAR! 00415 if (sc.isFailure()) return sc; 00416 } 00417 00419 HitPipeline::iterator current; 00421 HitPipeline::iterator previous; 00422 00423 TimeStamp GapStop; 00424 00425 bool GapFound=false; 00426 bool TimeTest=false; 00427 00428 // 00429 // _________hhhhhh__hh_h____hhhhh_H ________________hhhh_H_hhh_Hhhhhhhhhhh___ 00430 // Chunk Start Stop|GapStart Stop 00431 // 00432 00433 TimeStamp BeginTime; 00434 TimeStamp EndTime; 00435 // Gap must be found and there is no further hits from lower stage will get into this gap. 00436 do { 00437 debug()<<"Currently hit pipeline have SimHits #: "<<m_hitPipeline.size()<<endreq; 00438 BeginTime=m_hitPipeline.begin()->first; 00439 EndTime= m_hitPipeline.rbegin()->first; 00440 debug()<<"Pipeline begin time: "<<BeginTime<<endreq; 00441 debug()<<"Pipeline end time: "<<EndTime<<endreq; 00442 00443 current = m_hitPipeline.begin(); 00444 // reset them to the beginning of hits pipeline 00445 m_chunkStart = current->first; 00446 m_chunkStop = current->first; 00447 00448 ++current; // Start from the second one 00449 for (;current != m_hitPipeline.end(); 00450 ++current) { // try all SimHits in one SimHeader 00451 previous = current; 00452 --previous; // Always compare current one with previous one 00453 double dt_sec = (current->first - previous->first); // in second 00454 double dt_units = dt_sec * CLHEP::second; // with unit 00455 00456 verbose() << "Previous hit at " << previous->first 00457 << " has deltaT " << dt_units/CLHEP::nanosecond << " ns" << endreq; 00458 verbose() << "Current hit at " << current->first << endreq; 00459 00461 if (dt_units >= m_minTimeGap) { 00462 m_chunkStop =previous->first; // the previous one's time, moving it forward 00463 GapStop=m_chunkStop; 00464 GapStop.Add(m_minTimeGap/CLHEP::second); // where the gap stop 00465 00466 GapFound=true; 00467 debug() << "Found gap of " << dt_units/CLHEP::nanosecond<< " ns" << endreq; 00468 if (lowerStage()->currentTime() >= GapStop) { 00469 TimeTest=true; 00470 debug() << "TimeTest true" <<endreq; 00471 debug() << "chunk start " << m_chunkStart <<endreq; 00472 debug() << "chunk stop " << m_chunkStop <<endreq; 00473 return StatusCode::SUCCESS; 00474 } else { 00475 TimeTest=false; 00476 } 00477 break; // Leave the innermost for loop 00478 } 00479 } // All SimHits in pipeline 00480 // Try some hits even gap is found. By asking more hits from lower stage the time can be updated. 00481 if ( !GapFound || !TimeTest ) { 00482 StatusCode sc = this->moar(); // need MOAR! 00483 if (sc.isFailure()) return sc; 00484 debug()<<"Right now hit pipeline have SimHits #: "<<m_hitPipeline.size()<<endreq; 00485 GapFound=false; 00486 TimeTest=false; 00487 } 00488 } while( !GapFound || !TimeTest ); //Both need to be satisfied. 00489 00490 return StatusCode::FAILURE; // It will never get here 00491 }
StatusCode ElecSimProc::fillHitHeader | ( | ) | [private] |
Definition at line 493 of file ElecSimProc.cc.
00494 { 00495 // Take SimHitHeaders from m_currentHeaders and fill hit 00496 // collections of fake amalSimHitHeader with all hits that have times 00497 // before chunk stop. 00498 // A lot of bookkeeping 00499 debug()<<"Fill an amalgamated SimHitHeader"<<endreq; 00500 m_count=0; // count for how many hits add to the fake SimHitHeader 00501 m_simHeaderCount.clear(); 00502 00503 m_realHitEarliest = TimeStamp::GetEOT(); 00504 m_realHitLatest = TimeStamp::GetBOT(); 00505 00506 // Clear out any previously filled hit collections 00507 DayaBay::SimHitHeader::hc_map& amal_hcmap = m_amalSimHitHeader.hitCollection(); // detectors collection 00508 DayaBay::SimHitHeader::hc_map::const_iterator it, done = amal_hcmap.end(); 00509 for (it = amal_hcmap.begin(); it != done; ++it) { 00510 DayaBay::SimHitCollection* pHitCollection = it->second; 00511 pHitCollection->collection().clear(); 00512 delete pHitCollection; 00513 } 00514 amal_hcmap.clear(); 00515 00516 // Loop over each previously saved header that contributed at 00517 // least one hit that came within trunk window to the pipeline. 00518 for (size_t ihead = 0; ihead < m_currentHeaders.size(); ++ihead) { 00519 m_simHeaderCount[m_currentHeaders[ihead]]=0; 00520 00521 TimeStamp referenceTime = m_currentHeaders[ihead]->timeStamp(); 00522 const DayaBay::SimHitHeader* shhead = m_currentHeaders[ihead]->hits(); 00523 const DayaBay::SimHitHeader::hc_map& hcMap = shhead->hitCollection(); 00524 DayaBay::SimHitHeader::hc_map::const_iterator hcit, hcdone = hcMap.end(); 00525 00526 // loop over each detectors hit collection 00527 for (hcit = hcMap.begin(); hcit != hcdone; ++hcit) { 00528 00529 // Get source hit collection for this detector 00530 DayaBay::SimHitCollection* shcol = hcit->second; 00531 if (!shcol) { return StatusCode::FAILURE; } 00532 const DayaBay::SimHitCollection::hit_container& hc = shcol->collection(); 00533 00534 // copy requisite hits 00535 for (size_t ihit=0; ihit < hc.size(); ++ihit) { 00536 DayaBay::SimHit* hit = hc[ihit]; 00537 TimeStamp ts = referenceTime; 00538 if( hit->hitTime()/CLHEP::second > m_longest ) continue; 00539 ts.Add(hit->hitTime()/CLHEP::second); 00540 00541 //debug()<<"reference time "<< referenceTime<<endreq; 00542 //debug()<<"delta t "<<hit->hitTime()/CLHEP::second<<endreq; 00543 00544 if (ts >= m_chunkStart && ts <= m_chunkStop) { 00545 m_count++; 00546 m_simHeaderCount[m_currentHeaders[ihead]]++; 00547 00548 // 00549 if( ts<m_realHitEarliest ) m_realHitEarliest = ts; 00550 if( ts>m_realHitLatest ) m_realHitLatest = ts; 00551 00552 // Get/make our hit collection for this detector 00553 DayaBay::SimHitHeader::hc_map::iterator it_det=amal_hcmap.find(hcit->first); 00554 DayaBay::SimHitCollection* my_shcol; 00555 if (it_det == amal_hcmap.end()) { 00556 my_shcol = new DayaBay::SimHitCollection(); 00557 //....... what to set its SimHitHeader to? ......... 00558 my_shcol->setDetector(DayaBay::Detector(hcit->first)); 00559 amal_hcmap[hcit->first] = my_shcol; 00560 } else { 00561 my_shcol=it_det->second; 00562 } 00563 00564 if (!ihead) { 00565 // Set SimHitHeader to be the one in the first SimHeader. 00566 // Is this the right thing to do? 00567 my_shcol->setHeader(const_cast<DayaBay::SimHitHeader*>(shhead)); 00568 } 00569 00570 DayaBay::SimHitCollection::hit_container& my_hc = my_shcol->collection(); 00571 verbose()<<"hit with time "<<ts<<" added to amalgamated SimHitHeader"<<endreq; 00572 my_hc.push_back(hit); 00573 } 00574 00575 } // hits 00576 } // detectors 00577 } // headers 00578 00579 debug()<<m_count<<" hits added to the amalgamated SimHitHeader"<<endreq; 00580 debug()<<"Size of m_amalSimHitHeader.hitCollection "<<m_amalSimHitHeader.hitCollection().size()<<endreq; 00581 debug()<<"Size of m_simHeaderCount "<<m_simHeaderCount.size()<<endreq; 00582 00583 return StatusCode::SUCCESS; 00584 }
StatusCode ElecSimProc::findGiantHitInTime | ( | ) | [private] |
Definition at line 588 of file ElecSimProc.cc.
00589 { 00590 m_giantHitsInTime.clear(); 00591 HitPipeline::iterator ci, ci_end = m_giantHits.end(); 00592 00593 for( ci=m_giantHits.begin(); ci!=ci_end; ++ci ) { 00594 00595 TimeStamp ts = ci->first; 00596 DayaBay::SimHit* aHit = ci->second; 00597 00598 if (ts >= m_chunkStart && ts <= m_chunkStop) { 00599 m_giantHitsInTime.insert( pair<TimeStamp,DayaBay::SimHit*>(ts,aHit) ); 00600 } 00601 } 00602 return StatusCode::SUCCESS; 00603 }
StatusCode ElecSimProc::ES_execute | ( | ) | [private] |
Definition at line 607 of file ElecSimProc.cc.
00608 { 00609 DayaBay::ElecHeader* elecHeader = MakeHeaderObject(); 00610 00611 // Set input SimHeaders who has really contribute to the elecHeader 00612 // set context as the first simheader 00613 bool found; 00614 found=false; 00615 for (size_t ind=0; ind<m_currentHeaders.size(); ++ind) { 00616 if(m_simHeaderCount[m_currentHeaders[ind]]>0) { // at least 1 hit is in use. 00617 this->AppendInputHeader(m_currentHeaders[ind]); 00618 if(!found) { 00619 elecHeader->setContext(m_currentHeaders[ind]->context()); 00620 found=true; 00621 } 00622 } 00623 } 00624 00626 // TimeStamps 00627 debug()<<"Before setContext "<<elecHeader->timeStamp()<<endreq; 00628 00629 m_realHitEarliest.Add(-DayaBay::preTimeTolerance/CLHEP::second); // earliest possible time 00630 m_realHitLatest.Add(DayaBay::postTimeTolerance/CLHEP::second); // latest possible time 00631 00632 elecHeader->setEarliest(m_realHitEarliest); 00633 elecHeader->setLatest(m_realHitLatest); 00634 elecHeader->setTimeStamp(m_realHitEarliest); 00635 00636 m_realElecHeaderTime=m_realHitEarliest; // They should be always be the same. 00637 debug()<<"After setContext "<<elecHeader->timeStamp()<<endreq; 00638 00639 // Add the pulse collection header 00640 DayaBay::ElecPulseHeader* pulseHeader = 00641 new DayaBay::ElecPulseHeader(elecHeader); 00642 elecHeader->setPulseHeader(pulseHeader); 00643 00644 // Add the crate header 00645 DayaBay::ElecCrateHeader* crateHeader = 00646 new DayaBay::ElecCrateHeader(elecHeader); 00647 elecHeader->setCrateHeader(crateHeader); 00648 00649 // Send each detector's hit collections through ElecSim tool chains 00650 // It is from m_amalSimHitHeader. 00651 const DayaBay::SimHitHeader::hc_map& hcMap = m_amalSimHitHeader.hitCollection(); 00652 00653 debug() << "Processing hit collections" <<endreq; 00654 00656 DayaBay::SimHitHeader::hc_map::const_iterator hcIter; 00657 for (hcIter=hcMap.begin(); hcIter != hcMap.end(); ++hcIter) { 00658 00659 DayaBay::Detector det(hcIter->first); 00660 debug() << "Checking " << det.detName() << " for hits." << endreq; 00661 DayaBay::SimHitCollection* hits = hcIter->second; 00662 if(!hits) return StatusCode::FAILURE; 00663 00664 debug() << "Got hit collection from " << det.detName()<<" (id= " 00665 << det.siteDetPackedData() << ") " << " with " 00666 << hits->collection().size() << " hits." << endreq; 00667 00668 StatusCode sc; 00669 // Check if this detector should be simulated. 00670 if( std::find(m_detectors.begin(),m_detectors.end(),det) 00671 == m_detectors.end()){ 00672 debug() << "Detector " << det.detName() << " not simulated." << endreq; 00673 continue; 00674 } 00675 if( det.isAD() || det.isWaterShield() ){ 00676 // Process PMT hits 00677 debug() << "Processing PMT hits" <<endreq; 00678 DayaBay::ElecPulseCollection* pulses = 00679 new DayaBay::ElecPulseCollection(pulseHeader, det); 00680 sc = m_pmtTool->generatePulses(hits, pulses); 00681 if( sc != StatusCode::SUCCESS ) return sc; 00682 pulseHeader->addPulseCollection(pulses); 00683 00684 DayaBay::ElecFeeCrate* crate = new DayaBay::ElecFeeCrate(det, 00685 crateHeader); 00686 sc = m_feeTool->generateSignals(pulses, crate); 00687 if( sc != StatusCode::SUCCESS ) return sc; 00688 crateHeader->addCrate(crate); 00689 00690 }else if(det.detectorId() == DetectorId::kRPC){ 00691 // Process RPC hits 00692 debug() << "Processing RPC hits" <<endreq; 00693 DayaBay::ElecPulseCollection* pulses = 00694 new DayaBay::ElecPulseCollection(pulseHeader, det); 00695 sc = m_rpcTool->generatePulses(hits, pulses); 00696 if( sc != StatusCode::SUCCESS ) return sc; 00697 pulseHeader->addPulseCollection(pulses); 00698 00699 DayaBay::ElecFecCrate* crate = new DayaBay::ElecFecCrate(det, 00700 crateHeader); 00701 sc = m_fecTool->generateSignals(pulses, crate); 00702 if( sc != StatusCode::SUCCESS ) return sc; 00703 crateHeader->addCrate(crate); 00704 }else{ 00705 error() << "Unknown detector " << det << endreq; 00706 return StatusCode::FAILURE; 00707 } 00708 } 00709 00710 // Send result to stage. 00711 ElecData* ed = new ElecData(*elecHeader); 00712 thisStage()->pushElement(ed); 00713 this->registerData(*ed); 00714 debug()<<"to grep: (Full Elec) new data pushed out at time "<<m_realHitEarliest<<endreq; 00715 00716 return StatusCode::SUCCESS; 00717 }
StatusCode ElecSimProc::fastElecSim | ( | ) | [private] |
Definition at line 793 of file ElecSimProc.cc.
00794 { 00795 StatusCode sc; 00796 00797 TimeStamp TimeOfLast = TimeStamp::GetBOT(); 00799 HitPipeline::iterator it, it_end = m_giantHitsInTime.end(); 00800 00801 for( it=m_giantHitsInTime.begin(); it!=it_end; ++it ) { 00802 00803 sc = this->preExecute(); 00804 if (sc.isFailure()) return sc; 00805 00806 DayaBay::ElecHeader* elecHeader = MakeHeaderObject(); 00807 DayaBay::SimHit* aHit = it->second; 00808 TimeStamp ts = it->first; 00809 00810 TimeStamp earliest = ts; 00811 TimeStamp latest = ts; 00812 earliest.Add(-DayaBay::preTimeTolerance/CLHEP::second); // earliest possible time 00813 earliest.Add(-1e-9); // 1ns earlier than its real simulated products 00814 // Then it can always be popped out and simulated first 00815 latest.Add(DayaBay::postTimeTolerance/CLHEP::second); 00816 00817 //const DayaBay::SimHeader* sheader = m_giantHitSheader[aHit]; 00818 map<const DayaBay::SimHit*, const DayaBay::SimHeader*>::iterator itmap = m_giantHitSheader.find(aHit); 00819 const DayaBay::SimHeader* sheader; 00820 if( itmap != m_giantHitSheader.end() ) { 00821 sheader = itmap->second; 00822 } else { 00823 error()<<"Can't find aGiantHit's parent"<<endreq; 00824 } 00825 00826 debug()<<"earliest = "<<earliest<<endreq; 00827 debug()<<"aHit at "<<aHit<<" SimHeader at "<<sheader<<endreq; 00828 00830 elecHeader->setContext(sheader->context()); 00832 elecHeader->setEarliest(earliest); 00833 elecHeader->setLatest(latest); 00834 elecHeader->setTimeStamp(earliest); 00835 00836 if( earliest > TimeOfLast ) { 00837 TimeOfLast = earliest; 00838 } 00839 00841 this->AppendInputHeader( sheader ); 00842 00843 // Send result to stage. 00844 thisStage()->pushElement(new ElecData(*elecHeader)); 00845 debug()<<"to grep: (Fast Elec) new data pushed out at time "<<earliest<<endreq; 00846 debug()<<"elecHeader at "<<elecHeader<<endreq; 00847 // Add ElecPulseHeader and ElecCrateHeader 00848 DayaBay::ElecPulseHeader* pPulse = 0; // new DayaBay::ElecPulseHeader; 00849 DayaBay::ElecCrateHeader* pCrate = 0; // new DayaBay::ElecCrateHeader; 00850 elecHeader->setPulseHeader(pPulse); 00851 elecHeader->setCrateHeader(pCrate); 00852 00853 sc = this->postExecute(); 00854 if (sc.isFailure()) return sc; 00855 } 00856 00857 if( m_emptyElecHeaderTime < TimeOfLast ) { 00858 m_emptyElecHeaderTime = TimeOfLast; 00859 } else { 00860 error()<<"New simulated data has earlier time. Why?"<<endreq; 00861 } 00862 00863 return StatusCode::SUCCESS; 00864 }
StatusCode ElecSimProc::blowChunks | ( | ) | [private] |
Definition at line 719 of file ElecSimProc.cc.
00720 { 00721 // clean pipeline, remove all the hits before m_chunkStop 00722 HitPipeline::iterator start = m_hitPipeline.begin(); 00723 HitPipeline::iterator end = m_hitPipeline.end(); 00724 HitPipeline::iterator stop = m_hitPipeline.end(); // chunkStop 00725 HitPipeline::iterator it; 00726 for( it=start; it!=end; ++it ) { 00727 if(it->first>m_chunkStop) { 00728 stop=it; 00729 break; // Leave the innermost for loop 00730 } 00731 } 00732 //m_hitPipeline.find(m_chunkStop); // only find the first element in line, not good 00733 m_hitPipeline.erase(start,stop); // [start, stop) are erased and stop is not included. 00734 00736 start = m_giantHits.begin(); 00737 end = m_giantHits.end(); 00738 stop = m_giantHits.end(); // chunkStop 00739 00740 for( it=start; it!=end; ++it ) { 00741 00742 DayaBay::SimHit* aHit = it->second; 00743 TimeStamp ts = it->first; 00744 00745 if( ts>m_chunkStop ) { 00746 stop = it; 00747 break; 00748 } 00749 00750 m_giantHitSheader.erase( aHit ); 00751 delete aHit; 00752 } 00753 m_giantHits.erase(start,stop); // [start,stop) 00754 00755 00756 // clean used SimHeader if it has no hit with time less than m_chunkStop 00757 std::vector<const DayaBay::SimHeader*>::iterator ihd; 00758 for(ihd=m_currentHeaders.begin();ihd!=m_currentHeaders.end();++ihd) { 00759 00760 // get the real latest 00761 TimeStamp latest=TimeStamp::GetEOT(); 00762 map<const DayaBay::SimHeader*, TimeStamp>::iterator it = m_realLatest.find(*ihd); 00763 if( it!=m_realLatest.end() ) { 00764 latest = it->second; 00765 } else { 00766 error()<<"Can't find the latest time of a SimHeader object"<<endreq; 00767 } 00768 // use the real latest to remove header objects 00769 if( latest<=m_chunkStop ) { 00770 // release the reference to the SimHeader 00771 DayaBay::SimHeader* nonconst = const_cast<DayaBay::SimHeader*>(*ihd); 00772 nonconst->release(); 00773 // remove its local copy 00774 m_currentHeaders.erase(ihd); 00775 ihd--; 00776 m_realLatest.erase(it); 00777 } 00778 } 00779 00780 { 00781 debug()<<"Current m_giantHitSheader end"<<endreq; 00782 map<const DayaBay::SimHit*, const DayaBay::SimHeader*>::iterator it,itend = m_giantHitSheader.end(); 00783 for( it=m_giantHitSheader.begin(); it!=itend; it++ ) { 00784 debug()<<"aGiantHit "<<it->first<<" with SimHeader "<<it->second<<endreq; 00785 } 00786 } 00787 00788 return StatusCode::SUCCESS; 00789 }
StatusCode StageProcessor< DayaBay::ElecHeader >::registerData | ( | IStageData & | data | ) | [inherited] |
IStage * StageProcessor< DayaBay::ElecHeader >::thisStage | ( | ) | [inherited] |
IStage * StageProcessor< DayaBay::ElecHeader >::lowerStage | ( | ) | [inherited] |
std::vector<std::string> ElecSimProc::m_detectorNames [private] |
Definition at line 89 of file ElecSimProc.h.
std::string ElecSimProc::m_pmtToolName [private] |
Definition at line 91 of file ElecSimProc.h.
std::string ElecSimProc::m_rpcToolName [private] |
Definition at line 93 of file ElecSimProc.h.
std::string ElecSimProc::m_feeToolName [private] |
Definition at line 95 of file ElecSimProc.h.
std::string ElecSimProc::m_fecToolName [private] |
Definition at line 97 of file ElecSimProc.h.
std::vector<DayaBay::Detector> ElecSimProc::m_detectors [private] |
Definition at line 100 of file ElecSimProc.h.
IEsPulseTool* ElecSimProc::m_pmtTool [private] |
Definition at line 102 of file ElecSimProc.h.
IEsPulseTool* ElecSimProc::m_rpcTool [private] |
Definition at line 104 of file ElecSimProc.h.
IEsFrontEndTool* ElecSimProc::m_feeTool [private] |
Definition at line 106 of file ElecSimProc.h.
IEsFrontEndTool* ElecSimProc::m_fecTool [private] |
Definition at line 108 of file ElecSimProc.h.
FFTimeStamp ElecSimProc::m_currentTime [private] |
Definition at line 155 of file ElecSimProc.h.
int ElecSimProc::m_count [private] |
Definition at line 156 of file ElecSimProc.h.
std::vector<const DayaBay::SimHeader*> ElecSimProc::m_currentHeaders [private] |
Definition at line 159 of file ElecSimProc.h.
std::map<const DayaBay::SimHeader*,int> ElecSimProc::m_simHeaderCount [private] |
Definition at line 163 of file ElecSimProc.h.
std::map<const DayaBay::SimHeader*, TimeStamp > ElecSimProc::m_realLatest [private] |
Definition at line 169 of file ElecSimProc.h.
HitPipeline ElecSimProc::m_hitPipeline [private] |
Definition at line 174 of file ElecSimProc.h.
double ElecSimProc::m_minTimeGap [private] |
Definition at line 177 of file ElecSimProc.h.
TimeStamp ElecSimProc::m_chunkStart [private] |
Definition at line 181 of file ElecSimProc.h.
TimeStamp ElecSimProc::m_chunkStop [private] |
Definition at line 182 of file ElecSimProc.h.
TimeStamp ElecSimProc::m_realHitEarliest [private] |
Definition at line 184 of file ElecSimProc.h.
TimeStamp ElecSimProc::m_realHitLatest [private] |
Definition at line 185 of file ElecSimProc.h.
HitPipeline ElecSimProc::m_giantHits [private] |
A giant hit representing the geant4-noble partilce's time.
Definition at line 188 of file ElecSimProc.h.
std::map<const DayaBay::SimHit*, const DayaBay::SimHeader*> ElecSimProc::m_giantHitSheader [private] |
Keep track all parent SimHeader for created SimHit, just for bookkeeping.
Definition at line 190 of file ElecSimProc.h.
HitPipeline ElecSimProc::m_giantHitsInTime [private] |
TimeStamp ElecSimProc::m_realElecHeaderTime [private] |
The earliest time from real ElecHeader and empty ElecHeader.
Definition at line 195 of file ElecSimProc.h.
TimeStamp ElecSimProc::m_emptyElecHeaderTime [private] |
Definition at line 196 of file ElecSimProc.h.
double ElecSimProc::m_longest [private] |