/search.css" rel="stylesheet" type="text/css"/> /search.js">
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

In This Package:

Public Member Functions | Private Types | Private Member Functions | Private Attributes
ElecSimProc Class Reference

Fifteen stuff. More...

#include <ElecSimProc.h>

Inheritance diagram for ElecSimProc:
Inheritance graph
[legend]
Collaboration diagram for ElecSimProc:
Collaboration graph
[legend]

List of all members.

Public Member Functions

 ElecSimProc (const std::string &name, ISvcLocator *pSvcLocator)
virtual ~ElecSimProc ()
virtual StatusCode initialize ()
virtual StatusCode execute ()
StatusCode registerData (IStageData &data)
IStagethisStage ()
IStagelowerStage ()

Private Types

typedef HeaderStageData
< DayaBay::HeaderObject
HeaderData
typedef HeaderStageData
< DayaBay::SimHeader
SimData
typedef HeaderStageData
< DayaBay::ElecHeader
ElecData
typedef std::multimap
< TimeStamp, DayaBay::SimHit
*, std::less< TimeStamp > > 
HitPipeline

Private Member Functions

StatusCode getSimData (DayaBay::SimHeader *&output)
StatusCode fillPipeLine (DayaBay::SimHeader *shead)
StatusCode moar ()
StatusCode findGap ()
StatusCode fillHitHeader ()
StatusCode findGiantHitInTime ()
StatusCode ES_execute ()
StatusCode fastElecSim ()
StatusCode blowChunks ()

Private Attributes

std::vector< std::string > m_detectorNames
std::string m_pmtToolName
std::string m_rpcToolName
std::string m_feeToolName
std::string m_fecToolName
std::vector< DayaBay::Detectorm_detectors
IEsPulseToolm_pmtTool
IEsPulseToolm_rpcTool
IEsFrontEndToolm_feeTool
IEsFrontEndToolm_fecTool
FFTimeStamp m_currentTime
 The earliest (smallest) time which this has provided.
DayaBay::SimHitHeader m_amalSimHitHeader
int m_count
std::vector< const
DayaBay::SimHeader * > 
m_currentHeaders
std::map< const
DayaBay::SimHeader *, int > 
m_simHeaderCount
std::map< const
DayaBay::SimHeader
*, TimeStamp
m_realLatest
HitPipeline m_hitPipeline
double m_minTimeGap
TimeStamp m_chunkStart
TimeStamp m_chunkStop
TimeStamp m_realHitEarliest
TimeStamp m_realHitLatest
HitPipeline m_giantHits
 A giant hit representing the geant4-noble partilce's time.
std::map< const
DayaBay::SimHit *, const
DayaBay::SimHeader * > 
m_giantHitSheader
 Keep track all parent SimHeader for created SimHit, just for bookkeeping.
HitPipeline m_giantHitsInTime
 The giant hit in current time window.
TimeStamp m_realElecHeaderTime
 The earliest time from real ElecHeader and empty ElecHeader.
TimeStamp m_emptyElecHeaderTime
double m_longest
 Longest jump.

Detailed Description

Fifteen stuff.

A pull aware algorithm for gluing ElecSim tools into Staged processing.

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.


Member Typedef Documentation

Definition at line 116 of file ElecSimProc.h.

Definition at line 117 of file ElecSimProc.h.

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.


Constructor & Destructor Documentation

ElecSimProc::ElecSimProc ( const std::string &  name,
ISvcLocator *  pSvcLocator 
)

Definition at line 36 of file ElecSimProc.cc.

    : StageProcessor<DayaBay::ElecHeader>(name,pSvcLocator)
    , m_pmtTool(0)
    , m_rpcTool(0)
    , m_feeTool(0)
    , m_fecTool(0)
{

    m_pullMode=true;

    // ElecSim
    declareProperty("Detectors",m_detectorNames,"List of active detectors");
    declareProperty("PmtTool",m_pmtToolName="EsPmtEffectPulseTool",
                    "Name of the PMT simulation tool");
    declareProperty("RpcTool",m_rpcToolName="EsIdealPulseTool",
                    "Name of the RPC simulation tool");
    declareProperty("FeeTool",m_feeToolName="EsIdealFeeTool",
                    "Name of the PMT Front-end electronics simulation tool");
    declareProperty("FecTool",m_fecToolName="EsIdealFecTool",
                    "Name of the RPC Front-end electronics simulation tool");
    // Default configuration is DayaBay near site, ADs only.
    m_detectorNames.push_back("DayaBayAD1");
    m_detectorNames.push_back("DayaBayAD2");
    m_detectorNames.push_back("DayaBayIWS");
    m_detectorNames.push_back("DayaBayOWS");
    m_detectorNames.push_back("DayaBayRPC");
    m_detectorNames.push_back("LingAoAD1");
    m_detectorNames.push_back("LingAoAD2");
    m_detectorNames.push_back("LingAoIWS");
    m_detectorNames.push_back("LingAoOWS");
    m_detectorNames.push_back("LingAoRPC");
    m_detectorNames.push_back("FarAD1");
    m_detectorNames.push_back("FarAD2");
    m_detectorNames.push_back("FarAD3");
    m_detectorNames.push_back("FarAD4");
    m_detectorNames.push_back("FarIWS");
    m_detectorNames.push_back("FarOWS");
    m_detectorNames.push_back("FarRPC");

    // special for this package
    declareProperty("MinGapTime", m_minTimeGap =
                    DayaBay::preTimeTolerance + DayaBay::postTimeTolerance,
                    "Minimum width for a gap that defines a chunk of data.");
    declareProperty("LongestJump", m_longest = 10*365*24*60*60, 
                    "A protection for longest jump in time stamp. Default is set to 10 year.");

    m_hitPipeline.clear();
    m_giantHits.clear();
    m_giantHitsInTime.clear();
    m_giantHitSheader.clear();

    m_currentTime = TimeStamp::GetBOT();
    m_realElecHeaderTime = TimeStamp::GetBOT();
    m_emptyElecHeaderTime = TimeStamp::GetBOT();
}
ElecSimProc::~ElecSimProc ( ) [virtual]

Definition at line 92 of file ElecSimProc.cc.

{
}

Member Function Documentation

StatusCode ElecSimProc::initialize ( ) [virtual]

ElecSim initialization

Reimplemented from StageProcessor< DayaBay::ElecHeader >.

Definition at line 96 of file ElecSimProc.cc.

{
    StatusCode sc = this->StageProcessor<DayaBay::ElecHeader>::initialize();
    if (sc.isFailure()) return sc;

    // Convert detector names to Detector IDs
    int nDetNames = m_detectorNames.size();
    for(int detIdx=0; detIdx<nDetNames; detIdx++){
        DayaBay::Detector det(m_detectorNames[detIdx]);
        if(det.site() == Site::kUnknown
           || det.detectorId() == DetectorId::kUnknown){
            error() << "Invalid detector name: " << m_detectorNames[detIdx] << endreq;
            return StatusCode::FAILURE;
        }
        m_detectors.push_back(det);
    }    

    // Get PMT Pulse tool
    try {
        m_pmtTool = tool<IEsPulseTool>(m_pmtToolName);
    }
    catch(const GaudiException& exg) {
        fatal() << "Failed to get pmt tool: \"" << m_pmtToolName << "\"" << endreq;
        return StatusCode::FAILURE;
    }
    info () << "Added tool " << m_pmtToolName << endreq;
  
    // Get RPC Pulse tool
    if(m_pmtToolName == m_rpcToolName){
        m_rpcTool = m_pmtTool;
    }else {
        try {
            m_rpcTool = tool<IEsPulseTool>(m_rpcToolName);
        }
        catch(const GaudiException& exg) {
            fatal() << "Failed to get rpc tool: \"" << m_rpcToolName << "\"" << endreq;
            return StatusCode::FAILURE;
        }
        info () << "Added tool " << m_rpcToolName << endreq;
    }

    // Get PMT Front-end Electronics tool
    try {
        m_feeTool = tool<IEsFrontEndTool>(m_feeToolName);
    }
    catch(const GaudiException& exg) {
        fatal() << "Failed to get fee tool: \"" << m_feeToolName << "\"" << endreq;
        return StatusCode::FAILURE;
    }
    info () << "Added tool " << m_feeToolName << endreq;

    // Get RPC Front-end Electronics tool
    try {
        m_fecTool = tool<IEsFrontEndTool>(m_fecToolName);
    }
    catch(const GaudiException& exg) {
        fatal() << "Failed to get fec tool: \"" << m_fecToolName << "\"" << endreq;
        return StatusCode::FAILURE;
    }
    info () << "Added tool " << m_fecToolName << endreq;
    
    // Set the 'partition' for this run in RunDataSvc
    IRunDataSvc* runDataSvc = svc<IRunDataSvc>("RunDataSvc",true);
    if(!runDataSvc){
        error() << "Failed to get RunDataSvc" << endreq;
        return StatusCode::FAILURE;
    }
    int task = 1;  // Tells RunDataSvc to use new simulation run data.
    Context emptyContext; // Need a dummy context for this service call
    ServiceMode svcMode(emptyContext, task);
    const DayaBay::RunData* runData = runDataSvc->runData(svcMode);
    if(!runData){
        error() << "Failed to get RunData for simulation" << endreq;
        return StatusCode::FAILURE;
    }
    DayaBay::RunData modifiedRunData(*runData);
    modifiedRunData.setDetectors(m_detectors);
    sc = runDataSvc->setRunData(modifiedRunData);
    if( !sc.isSuccess() ){
        error() << "Failed to set detectors in run" << endreq;
    }

    info()<<"minTimeGap "<<m_minTimeGap << " ns" <<endreq;

    return sc;
}
StatusCode ElecSimProc::execute ( ) [virtual]

Definition at line 184 of file ElecSimProc.cc.

{
    debug()<<"Executing ..."<<endreq;
    // Only do anything if the current stage time has advanced beyond
    // the last time we ran.
    static bool first_time = true;
    if (!first_time && m_currentTime > thisStage()->currentTime()) {
        return StatusCode::SUCCESS;
    }
    first_time = false;

    StatusCode sc;

    sc = this->findGap();
    if (sc.isFailure()) return sc;

    sc = this->fillHitHeader();
    if (sc.isFailure()) return sc;    

    sc = this->findGiantHitInTime();
    if (sc.isFailure()) return sc;
    
    // Number of normal hits is above zero
    if( m_count>0 )  {
        sc = this->preExecute();
        if (sc.isFailure()) return sc;
        sc = this->ES_execute();
        if (sc.isFailure()) return sc;
        sc = this->postExecute();
        if (sc.isFailure()) return sc;
    }
    
    // Number of giant hits is above zero
    if( m_giantHitsInTime.size()>0 )  {
        sc = this->fastElecSim();
        if (sc.isFailure()) return sc;
    }

    // update m_currentTime
    if( m_realElecHeaderTime > m_emptyElecHeaderTime )  {
        m_currentTime = m_realElecHeaderTime;
    } else {
        m_currentTime = m_emptyElecHeaderTime;
    }

    // blowChunks() must be the last one, because it will hold the last refcount of some
    // input headers. Before do this, it must be added to inputheaders in postExecute().
    sc = this->blowChunks();
    if (sc.isFailure()) return sc;    

    return StatusCode::SUCCESS;
}
StatusCode ElecSimProc::getSimData ( DayaBay::SimHeader *&  output) [private]

Try SimData first, if it is from DetSimProc

Or it is from LoadingProc

Definition at line 238 of file ElecSimProc.cc.

{
    // Currently we get it from the lower stage, but this function
    // will eventually handle getting it from AES when the SimHeaders
    // are read in from file.

    IStageData* stageData = 0;
    StatusCode sc = lowerStage()->nextElement(stageData);
    if (sc.isFailure()) {
        error() << "Failed to get data from detector simulation stage"
                << endreq;
        return sc;
    }
    
    DayaBay::SimHeader* shead = 0;
    
    SimData* input=0;
    input = dynamic_cast<SimData*>(stageData);
    
    if( input!=0 ) {
        shead = &input->header();
    } else {
        HeaderData* pHeaderData = dynamic_cast<HeaderData*>(stageData);
        DayaBay::HeaderObject* pHeaderObj = &pHeaderData->header();
        shead = dynamic_cast<DayaBay::SimHeader*>(pHeaderObj);
    }
    
    if( shead==0 ) {
        error() << "Failed to get SimHeader from lower stage" <<endreq;
    }

    debug() << "stage data's time "<< stageData->time()<<endreq;

    // we must temporarily hold on to this DataObject because this
    // input may have been created during previous exection cycle.
    shead->addRef();            // And because deleting the SimData
    delete stageData;           // drops the only other reference.  It
    stageData=0;                // will be released once it is no
                                // longer providing any hits to the
                                // pipeline

    debug() << "Current time of lower stage "<<lowerStage()->currentTime()<<endreq;
    debug() << "Got SimHeader from lower stage at " << shead->timeStamp() << endreq;

    output = shead;
    return StatusCode::SUCCESS;
}
StatusCode ElecSimProc::fillPipeLine ( DayaBay::SimHeader shead) [private]

Loop over all SimHitCollections

Loop over all GenHeader particles See if there is any interesting geant4-noble particles If there is, make a giant hit for them. (i.e. MuonProphet muon) One SimHeader, one GenHeader

Create a GiantHit for all geant4-noble particles

To do: set time based on the closest approach from track to detector right now, 0;

use weight 0 and kUnknown detector to identified them

Before go into real electronic simulation, this will be taken out.

Definition at line 290 of file ElecSimProc.cc.

{
    // pipeline
    TimeStamp realLatest = TimeStamp::GetBOT();

    // Get the SimHit Header
    const DayaBay::SimHitHeader* shitHeader = shead->hits();
    if (!shitHeader) {
        warning() << "Got NULL SimHitHeader with SimHeader @ "
                  << shead->timeStamp()
                  << " skipping..."
                  << endreq;
        return StatusCode::SUCCESS;
    }
    const DayaBay::SimHitHeader::hc_map& hcMap = shitHeader->hitCollection();

    const TimeStamp& absTime = shead->timeStamp();
    debug()<<"absTime = "<<absTime<<endreq;
    DayaBay::SimHitHeader::hc_map::const_iterator hcIter;
    for (hcIter=hcMap.begin(); hcIter != hcMap.end(); ++hcIter) {

        DayaBay::Detector det(hcIter->first);
        DayaBay::SimHitCollection* hits = hcIter->second;
        if(!hits) return StatusCode::FAILURE;
        
        const DayaBay::SimHitCollection::hit_container& hc = hits->collection();
        for (size_t ind=0; ind < hc.size(); ++ind) {
            DayaBay::SimHit* hit = hc[ind];

            TimeStamp ts = absTime;
            if( hit->hitTime()/CLHEP::second < m_longest ) {
                ts.Add(hit->hitTime()/CLHEP::second);
                m_hitPipeline.insert(pair<TimeStamp,DayaBay::SimHit*>(ts,hit));
                //debug()<<hit->hitTime()/CLHEP::second<<" [sec]"<<endreq;
                // find the real latest time
                if( ts>realLatest ) realLatest = ts;
            } else {
                warning()<<"A hit in the far future "<<hit->hitTime()/CLHEP::second<<" [sec] is skipped."<<endreq;
            }
        }
    }

    const DayaBay::IHeader* iheader = (shead->inputHeaders()) [0]; 
    const DayaBay::GenHeader* genHeader = dynamic_cast<const DayaBay::GenHeader*>(iheader);
    const HepMC::GenEvent* genEvent = genHeader->event();

    if( genEvent )  {
        HepMC::GenEvent::particle_const_iterator pci, pci_end = genEvent->particles_end();
    
        for( pci = genEvent->particles_begin(); pci != pci_end; ++pci )  {
            if( (*pci)->pdg_id() == 13 || (*pci)->pdg_id() == -13 )  {   // muon
                if( (*pci)->status() != MpMuonFate::kNeedSim )   {   //  geant4-noble particle
                    debug()<<"Find one geant4-noble particle"<<endreq;
                    DayaBay::SimHit* aGiantHit = new DayaBay::SimHit;
                    double dt = 0 * CLHEP::nanosecond;
                    aGiantHit->setHitTime(dt);
                    aGiantHit->setWeight(0);
                    aGiantHit->setSensDetId( DetectorId::kUnknown );
                    debug()<<"absTime "<<absTime<<endreq;
                    TimeStamp ts = absTime;
                    ts.Add(dt/CLHEP::second);
                    debug()<<"ts      "<<ts<<endreq;
                    m_hitPipeline.insert(pair<TimeStamp,DayaBay::SimHit*>(ts,aGiantHit));
                    m_giantHits.insert(pair<TimeStamp,DayaBay::SimHit*>(ts,aGiantHit));
                    m_giantHitSheader[aGiantHit] = shead;
                    debug()<<"aGiantHit is added "<<aGiantHit<<" with SimHeader "<<shead<<" at time "<<ts<<endreq;
                    debug()<<"Now m_giantHitSheader has "<<m_giantHitSheader.size()<<" giant hits."<<endreq;
                    // find the real latest time
                    if( ts>realLatest ) realLatest = ts;
                    break; // only process one fast simulated muon
                }
            }
        }
    }
    
    debug()<<"Current m_giantHitSheader begin"<<endreq;
    map<const DayaBay::SimHit*, const DayaBay::SimHeader*>::iterator it,itend = m_giantHitSheader.end();
    for( it=m_giantHitSheader.begin(); it!=itend; it++ ) {
        debug()<<"aGiantHit "<<it->first<<" with SimHeader "<<it->second<<endreq;
    }

    // current SimHeaders
    m_currentHeaders.push_back(shead);
    debug()<<"m_currentHeaders has "<<m_currentHeaders.size()<<" SimHeaders."<<endreq;

    // record the real latest time
    m_realLatest.insert( pair<DayaBay::SimHeader*, TimeStamp>(shead,realLatest) );
    
    return StatusCode::SUCCESS;
}
StatusCode ElecSimProc::moar ( ) [private]

Definition at line 391 of file ElecSimProc.cc.

{
    debug() << "Getting MOAR!" << endreq;

    DayaBay::SimHeader* shead = 0;
    StatusCode sc = this->getSimData(shead);
    if (sc.isFailure()) return sc;
    
    sc = this->fillPipeLine(shead);
    if (sc.isFailure()) return sc;

    return StatusCode::SUCCESS;
}
StatusCode ElecSimProc::findGap ( ) [private]

A iterator pointing to current one

A iterator pointint to the previous one

Gap found

Definition at line 405 of file ElecSimProc.cc.

{
    debug()<<"Entering findGap()"<<endreq;
    while ( m_hitPipeline.size() < 2 )  {  // at least two hits are needed to find a gap
        StatusCode sc = this->moar(); // need MOAR!
        if (sc.isFailure()) return sc;
    }

    HitPipeline::iterator current;
    HitPipeline::iterator previous;

    TimeStamp GapStop;

    bool GapFound=false;
    bool TimeTest=false;

    //
    // _________hhhhhh__hh_h____hhhhh_H ________________hhhh_H_hhh_Hhhhhhhhhhh___
    //          Chunk Start        Stop|GapStart    Stop
    //

    TimeStamp BeginTime;
    TimeStamp EndTime;
    // Gap must be found and there is no further hits from lower stage will get into this gap.
    do {
        debug()<<"Currently hit pipeline have SimHits #: "<<m_hitPipeline.size()<<endreq;
        BeginTime=m_hitPipeline.begin()->first;
        EndTime= m_hitPipeline.rbegin()->first;
        debug()<<"Pipeline begin time: "<<BeginTime<<endreq;
        debug()<<"Pipeline end time:   "<<EndTime<<endreq;

        current = m_hitPipeline.begin();
        // reset them to the beginning of hits pipeline
        m_chunkStart = current->first;
        m_chunkStop  = current->first;

        ++current; // Start from the second one
        for (;current != m_hitPipeline.end(); 
             ++current) {  // try all SimHits in one SimHeader
            previous = current;
            --previous;    // Always compare current one with previous one
            double dt_sec = (current->first - previous->first);  // in second
            double dt_units = dt_sec * CLHEP::second;            // with unit

            verbose() << "Previous hit at " << previous->first 
                      << " has deltaT " << dt_units/CLHEP::nanosecond << " ns" << endreq;
            verbose() << "Current hit at " << current->first << endreq;

            if (dt_units >= m_minTimeGap) {
                m_chunkStop =previous->first;             // the previous one's time, moving it forward
                GapStop=m_chunkStop;
                GapStop.Add(m_minTimeGap/CLHEP::second);  // where the gap stop

                GapFound=true;
                debug() << "Found gap of " << dt_units/CLHEP::nanosecond<< " ns" << endreq;
                if (lowerStage()->currentTime() >= GapStop) {
                    TimeTest=true;
                    debug() << "TimeTest true" <<endreq;
                    debug() << "chunk start " << m_chunkStart <<endreq;
                    debug() << "chunk stop  " << m_chunkStop <<endreq;
                    return StatusCode::SUCCESS;
                } else {
                    TimeTest=false;
                }               
                break;  // Leave the innermost for loop
            }
        } // All SimHits in pipeline
        // Try some hits even gap is found. By asking more hits from lower stage the time can be updated.
        if ( !GapFound || !TimeTest ) {
            StatusCode sc = this->moar(); // need MOAR!  
            if (sc.isFailure()) return sc;
            debug()<<"Right now hit pipeline have SimHits #: "<<m_hitPipeline.size()<<endreq;
            GapFound=false;
            TimeTest=false;
        }
    } while( !GapFound || !TimeTest ); //Both need to be satisfied.
     
    return StatusCode::FAILURE;  // It will never get here
}
StatusCode ElecSimProc::fillHitHeader ( ) [private]

Definition at line 488 of file ElecSimProc.cc.

{
    // Take SimHitHeaders from m_currentHeaders and fill hit
    // collections of fake amalSimHitHeader with all hits that have times
    // before chunk stop.
    // A lot of bookkeeping
    debug()<<"Fill an amalgamated SimHitHeader"<<endreq;
    m_count=0;  // count for how many hits add to the fake SimHitHeader
    m_simHeaderCount.clear();

    m_realHitEarliest = TimeStamp::GetEOT();
    m_realHitLatest   = TimeStamp::GetBOT();

    // Clear out any previously filled hit collections
    DayaBay::SimHitHeader::hc_map& amal_hcmap = m_amalSimHitHeader.hitCollection(); // detectors collection
    DayaBay::SimHitHeader::hc_map::const_iterator it, done = amal_hcmap.end();
    for (it = amal_hcmap.begin(); it != done; ++it) {
        DayaBay::SimHitCollection* pHitCollection = it->second;
        pHitCollection->collection().clear();
        delete pHitCollection;
    }
    amal_hcmap.clear();

    // Loop over each previously saved header that contributed at
    // least one hit that came within trunk window to the pipeline.
    for (size_t ihead = 0; ihead < m_currentHeaders.size(); ++ihead) {
        m_simHeaderCount[m_currentHeaders[ihead]]=0;
        
        TimeStamp referenceTime = m_currentHeaders[ihead]->timeStamp();
        const DayaBay::SimHitHeader* shhead = m_currentHeaders[ihead]->hits();
        const DayaBay::SimHitHeader::hc_map& hcMap = shhead->hitCollection();
        DayaBay::SimHitHeader::hc_map::const_iterator hcit, hcdone = hcMap.end();

        // loop over each detectors hit collection
        for (hcit = hcMap.begin(); hcit != hcdone; ++hcit) {

            // Get source hit collection for this detector
            DayaBay::SimHitCollection* shcol = hcit->second;
            if (!shcol) { return StatusCode::FAILURE; }
            const DayaBay::SimHitCollection::hit_container& hc = shcol->collection();
            
            // copy requisite hits
            for (size_t ihit=0; ihit < hc.size(); ++ihit) {
                DayaBay::SimHit* hit = hc[ihit];
                TimeStamp ts = referenceTime;
                if( hit->hitTime()/CLHEP::second > m_longest ) continue;
                ts.Add(hit->hitTime()/CLHEP::second);
                
                //debug()<<"reference time "<< referenceTime<<endreq;
                //debug()<<"delta t "<<hit->hitTime()/CLHEP::second<<endreq;

                if (ts >= m_chunkStart && ts <= m_chunkStop) {
                    m_count++;
                    m_simHeaderCount[m_currentHeaders[ihead]]++;

                    // 
                    if( ts<m_realHitEarliest ) m_realHitEarliest = ts;
                    if( ts>m_realHitLatest )   m_realHitLatest = ts;

                    // Get/make our hit collection for this detector
                    DayaBay::SimHitHeader::hc_map::iterator it_det=amal_hcmap.find(hcit->first);
                    DayaBay::SimHitCollection* my_shcol;
                    if (it_det == amal_hcmap.end()) {
                        my_shcol = new DayaBay::SimHitCollection();
                        //....... what to set its SimHitHeader to? ......... 
                        my_shcol->setDetector(DayaBay::Detector(hcit->first));
                        amal_hcmap[hcit->first] = my_shcol;
                    } else {
                        my_shcol=it_det->second;
                    }
                    
                    if (!ihead) {
                        // Set SimHitHeader to be the one in the first SimHeader.
                        // Is this the right thing to do?
                        my_shcol->setHeader(const_cast<DayaBay::SimHitHeader*>(shhead));
                    }
                    
                    DayaBay::SimHitCollection::hit_container& my_hc = my_shcol->collection();
                    verbose()<<"hit with time "<<ts<<" added to amalgamated SimHitHeader"<<endreq;
                    my_hc.push_back(hit);
                }

            } // hits
        } // detectors
    } // headers

    debug()<<m_count<<" hits added to the amalgamated SimHitHeader"<<endreq;
    debug()<<"Size of m_amalSimHitHeader.hitCollection "<<m_amalSimHitHeader.hitCollection().size()<<endreq;
    debug()<<"Size of m_simHeaderCount "<<m_simHeaderCount.size()<<endreq;

    return StatusCode::SUCCESS;
}
StatusCode ElecSimProc::findGiantHitInTime ( ) [private]

Definition at line 583 of file ElecSimProc.cc.

{
    m_giantHitsInTime.clear();
    HitPipeline::iterator ci, ci_end = m_giantHits.end();
    
    for( ci=m_giantHits.begin(); ci!=ci_end; ++ci ) {
        
        TimeStamp ts = ci->first;
        DayaBay::SimHit* aHit = ci->second;
        
        if (ts >= m_chunkStart && ts <= m_chunkStop) {
            m_giantHitsInTime.insert( pair<TimeStamp,DayaBay::SimHit*>(ts,aHit) );
        }
    }
    return StatusCode::SUCCESS;
}
StatusCode ElecSimProc::ES_execute ( ) [private]

setTimestamp must be after setContext, otherwise it will be overwritten.

Loop over all SimHitCollections

Definition at line 602 of file ElecSimProc.cc.

{
    DayaBay::ElecHeader* elecHeader = MakeHeaderObject();

    // Set input SimHeaders who has really contribute to the elecHeader
    // set context as the first simheader
    bool found;
    found=false;
    for (size_t ind=0; ind<m_currentHeaders.size(); ++ind) {
        if(m_simHeaderCount[m_currentHeaders[ind]]>0) {  // at least 1 hit is in use.
            this->AppendInputHeader(m_currentHeaders[ind]);
            if(!found) {
                elecHeader->setContext(m_currentHeaders[ind]->context());
                found=true;
            }
        }
    }

    // TimeStamps
    debug()<<"Before setContext "<<elecHeader->timeStamp()<<endreq;

    m_realHitEarliest.Add(-DayaBay::preTimeTolerance/CLHEP::second);  // earliest possible time
    m_realHitLatest.Add(DayaBay::postTimeTolerance/CLHEP::second);    // latest possible time

    elecHeader->setEarliest(m_realHitEarliest);
    elecHeader->setLatest(m_realHitLatest);
    elecHeader->setTimeStamp(m_realHitEarliest);

    m_realElecHeaderTime=m_realHitEarliest;               // They should be always be the same.
    debug()<<"After setContext  "<<elecHeader->timeStamp()<<endreq;

    // Add the pulse collection header
    DayaBay::ElecPulseHeader* pulseHeader = 
        new DayaBay::ElecPulseHeader(elecHeader);
    elecHeader->setPulseHeader(pulseHeader);

    // Add the crate header
    DayaBay::ElecCrateHeader* crateHeader = 
        new DayaBay::ElecCrateHeader(elecHeader);
    elecHeader->setCrateHeader(crateHeader);

    // Send each detector's hit collections through ElecSim tool chains
    // It is from m_amalSimHitHeader.
    const DayaBay::SimHitHeader::hc_map& hcMap = m_amalSimHitHeader.hitCollection();

    debug() << "Processing hit collections" <<endreq;

    DayaBay::SimHitHeader::hc_map::const_iterator hcIter;
    for (hcIter=hcMap.begin(); hcIter != hcMap.end(); ++hcIter) {

        DayaBay::Detector det(hcIter->first);
        debug() << "Checking " << det.detName() << " for hits." << endreq;
        DayaBay::SimHitCollection* hits = hcIter->second;
        if(!hits) return StatusCode::FAILURE;

        debug() << "Got hit collection from " << det.detName()<<" (id= "
                << det.siteDetPackedData() << ") " << " with "
                << hits->collection().size() << " hits." << endreq;

        StatusCode sc;
        // Check if this detector should be simulated.
        if( std::find(m_detectors.begin(),m_detectors.end(),det)
            == m_detectors.end()){
            debug() << "Detector " << det.detName() << " not simulated." << endreq;
            continue;
        }
        if( det.isAD() || det.isWaterShield() ){
            // Process PMT hits
            debug() << "Processing PMT hits" <<endreq;
            DayaBay::ElecPulseCollection* pulses =
                new DayaBay::ElecPulseCollection(pulseHeader, det);
            sc = m_pmtTool->generatePulses(hits, pulses);
            if( sc != StatusCode::SUCCESS ) return sc;
            pulseHeader->addPulseCollection(pulses);

            DayaBay::ElecFeeCrate* crate = new DayaBay::ElecFeeCrate(det,
                                                                     crateHeader);
            sc = m_feeTool->generateSignals(pulses, crate);
            if( sc != StatusCode::SUCCESS ) return sc;
            crateHeader->addCrate(crate);

        }else if(det.detectorId() == DetectorId::kRPC){
            // Process RPC hits
            debug() << "Processing RPC hits" <<endreq;
            DayaBay::ElecPulseCollection* pulses =
                new DayaBay::ElecPulseCollection(pulseHeader, det);
            sc = m_rpcTool->generatePulses(hits, pulses);
            if( sc != StatusCode::SUCCESS ) return sc;
            pulseHeader->addPulseCollection(pulses);

            DayaBay::ElecFecCrate* crate = new DayaBay::ElecFecCrate(det,
                                                                     crateHeader);
            sc = m_fecTool->generateSignals(pulses, crate);
            if( sc != StatusCode::SUCCESS ) return sc;
            crateHeader->addCrate(crate);
        }else{
            error() << "Unknown detector " << det << endreq;
            return StatusCode::FAILURE;
        }
    }

    // Send result to stage.
    ElecData* ed = new ElecData(*elecHeader);
    thisStage()->pushElement(ed);
    this->registerData(*ed);
    debug()<<"to grep: (Full Elec) new data pushed out at time "<<m_realHitEarliest<<endreq;
    
    return StatusCode::SUCCESS;
}
StatusCode ElecSimProc::fastElecSim ( ) [private]

Each giant hit will have one empty ElecHeader

set context

set time

set inputHeaders

Definition at line 788 of file ElecSimProc.cc.

{
    StatusCode sc;

    TimeStamp TimeOfLast = TimeStamp::GetBOT();
    HitPipeline::iterator it, it_end = m_giantHitsInTime.end();

    for( it=m_giantHitsInTime.begin(); it!=it_end; ++it )  {
        
        sc = this->preExecute();
        if (sc.isFailure()) return sc;

        DayaBay::ElecHeader* elecHeader = MakeHeaderObject();
        DayaBay::SimHit* aHit = it->second;
        TimeStamp ts = it->first;

        TimeStamp earliest = ts;
        TimeStamp latest =   ts;
        earliest.Add(-DayaBay::preTimeTolerance/CLHEP::second);  // earliest possible time
        earliest.Add(-1e-9); // 1ns earlier than its real simulated products
                             // Then it can always be popped out and simulated first
        latest.Add(DayaBay::postTimeTolerance/CLHEP::second);
        
        //const DayaBay::SimHeader* sheader = m_giantHitSheader[aHit];
        map<const DayaBay::SimHit*, const DayaBay::SimHeader*>::iterator itmap = m_giantHitSheader.find(aHit);
        const DayaBay::SimHeader* sheader;
        if( itmap != m_giantHitSheader.end() ) {
            sheader = itmap->second;
        } else {
            error()<<"Can't find aGiantHit's parent"<<endreq;
        }

        debug()<<"earliest = "<<earliest<<endreq;
        debug()<<"aHit at "<<aHit<<"   SimHeader at "<<sheader<<endreq;

        elecHeader->setContext(sheader->context());
        elecHeader->setEarliest(earliest);
        elecHeader->setLatest(latest);
        elecHeader->setTimeStamp(earliest);
        
        if( earliest > TimeOfLast ) {
            TimeOfLast = earliest;
        }

        this->AppendInputHeader( sheader );

        // Send result to stage.
        thisStage()->pushElement(new ElecData(*elecHeader));
        debug()<<"to grep: (Fast Elec) new data pushed out at time "<<earliest<<endreq;
        debug()<<"elecHeader at "<<elecHeader<<endreq;
        // Add ElecPulseHeader and ElecCrateHeader
        DayaBay::ElecPulseHeader* pPulse = 0; // new DayaBay::ElecPulseHeader;
        DayaBay::ElecCrateHeader* pCrate = 0; // new DayaBay::ElecCrateHeader;
        elecHeader->setPulseHeader(pPulse);
        elecHeader->setCrateHeader(pCrate);

        sc = this->postExecute();
        if (sc.isFailure()) return sc;
    }

    if( m_emptyElecHeaderTime < TimeOfLast ) {
        m_emptyElecHeaderTime = TimeOfLast;
    } else {
        error()<<"New simulated data has earlier time. Why?"<<endreq;
    }
    
    return StatusCode::SUCCESS;
}
StatusCode ElecSimProc::blowChunks ( ) [private]

clean m_giantHits

Definition at line 714 of file ElecSimProc.cc.

{
    // clean pipeline, remove all the hits before m_chunkStop
    HitPipeline::iterator start = m_hitPipeline.begin();
    HitPipeline::iterator end   = m_hitPipeline.end();
    HitPipeline::iterator stop  = m_hitPipeline.end();  // chunkStop
    HitPipeline::iterator it;
    for( it=start; it!=end; ++it )  {
        if(it->first>m_chunkStop) {
            stop=it;
            break;  // Leave the innermost for loop
        }
    }
    //m_hitPipeline.find(m_chunkStop); // only find the first element in line, not good
    m_hitPipeline.erase(start,stop);    // [start, stop) are erased and stop is not included.

    start = m_giantHits.begin();
    end   = m_giantHits.end();
    stop  = m_giantHits.end(); // chunkStop

    for( it=start; it!=end; ++it )  {
        
        DayaBay::SimHit* aHit = it->second;
        TimeStamp ts = it->first;
        
        if( ts>m_chunkStop )  {
            stop = it;
            break;
        }
        
        m_giantHitSheader.erase( aHit );
        delete aHit;
    }
    m_giantHits.erase(start,stop);  // [start,stop)


    // clean used SimHeader if it has no hit with time less than m_chunkStop
    std::vector<const DayaBay::SimHeader*>::iterator ihd;
    for(ihd=m_currentHeaders.begin();ihd!=m_currentHeaders.end();++ihd) {

        // get the real latest
        TimeStamp latest=TimeStamp::GetEOT();
        map<const DayaBay::SimHeader*, TimeStamp>::iterator it = m_realLatest.find(*ihd);
        if( it!=m_realLatest.end() ) {
            latest = it->second;
        } else {
            error()<<"Can't find the latest time of a SimHeader object"<<endreq;
        }
        // use the real latest to remove header objects
        if( latest<=m_chunkStop ) {
            // release the reference to the SimHeader
            DayaBay::SimHeader* nonconst = const_cast<DayaBay::SimHeader*>(*ihd);
            nonconst->release();
            // remove its local copy
            m_currentHeaders.erase(ihd);
            ihd--;
            m_realLatest.erase(it);
        }
    }

    {
        debug()<<"Current m_giantHitSheader end"<<endreq;
        map<const DayaBay::SimHit*, const DayaBay::SimHeader*>::iterator it,itend = m_giantHitSheader.end();
        for( it=m_giantHitSheader.begin(); it!=itend; it++ ) {
            debug()<<"aGiantHit "<<it->first<<" with SimHeader "<<it->second<<endreq;
        }
    }

    return StatusCode::SUCCESS;
}
StatusCode StageProcessor< class >::registerData ( IStageData data) [inherited]
IStage * StageProcessor< class >::thisStage ( ) [inherited]
IStage * StageProcessor< class >::lowerStage ( ) [inherited]

Member Data Documentation

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.

Definition at line 100 of file ElecSimProc.h.

Definition at line 102 of file ElecSimProc.h.

Definition at line 104 of file ElecSimProc.h.

Definition at line 106 of file ElecSimProc.h.

Definition at line 108 of file ElecSimProc.h.

The earliest (smallest) time which this has provided.

Definition at line 111 of file ElecSimProc.h.

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.

Definition at line 169 of file ElecSimProc.h.

Definition at line 174 of file ElecSimProc.h.

double ElecSimProc::m_minTimeGap [private]

Definition at line 177 of file ElecSimProc.h.

Definition at line 181 of file ElecSimProc.h.

Definition at line 182 of file ElecSimProc.h.

Definition at line 184 of file ElecSimProc.h.

Definition at line 185 of file ElecSimProc.h.

A giant hit representing the geant4-noble partilce's time.

Definition at line 188 of file ElecSimProc.h.

Keep track all parent SimHeader for created SimHit, just for bookkeeping.

Definition at line 190 of file ElecSimProc.h.

The giant hit in current time window.

Definition at line 192 of file ElecSimProc.h.

The earliest time from real ElecHeader and empty ElecHeader.

Definition at line 195 of file ElecSimProc.h.

Definition at line 196 of file ElecSimProc.h.

double ElecSimProc::m_longest [private]

Longest jump.

Definition at line 199 of file ElecSimProc.h.


The documentation for this class was generated from the following files:
| Classes | Job Modules | Data Objects | Services | Algorithms | Tools | Packages | Directories | Tracs |

Generated on Fri May 16 2014 10:22:01 for ElecSimProc by doxygen 1.7.4