#include <CaloSD.h>

Inheritance diagram for CaloSD:

Public Member Functions
	CaloSD (G4String aSDname, const DDCompactView &cpv, const SensitiveDetectorCatalog &clg, edm::ParameterSet const &p, const SimTrackManager *, int tSlice=1, bool ignoreTkID=false)

virtual void	clear ()

virtual void	DrawAll ()

virtual void	EndOfEvent (G4HCofThisEvent *eventHC)

void	fillHits (edm::PCaloHitContainer &, std::string n)

virtual double	getEnergyDeposit (G4Step *step)

virtual void	Initialize (G4HCofThisEvent *HCE)

virtual void	PrintAll ()

virtual bool	ProcessHits (G4Step step, G4TouchableHistory tHistory)

virtual bool	ProcessHits (G4GFlashSpot aSpot, G4TouchableHistory )

virtual uint32_t	setDetUnitId (G4Step *step)=0

virtual	~CaloSD ()

Public Member Functions inherited from SensitiveCaloDetector
	SensitiveCaloDetector (std::string &iname, const DDCompactView &cpv, const SensitiveDetectorCatalog &clg, edm::ParameterSet const &p)

Public Member Functions inherited from SensitiveDetector
virtual void	AssignSD (const std::string &vname)

Local3DPoint	ConvertToLocal3DPoint (const G4ThreeVector &point)

Local3DPoint	FinalStepPosition (G4Step *s, coordinates)

virtual std::vector< std::string >	getNames ()

Local3DPoint	InitialStepPosition (G4Step *s, coordinates)

std::string	nameOfSD ()

void	NaNTrap (G4Step *step)

void	Register ()

	SensitiveDetector (std::string &iname, const DDCompactView &cpv, const SensitiveDetectorCatalog &, edm::ParameterSet const &p)

virtual	~SensitiveDetector ()

Public Member Functions inherited from Observer< const BeginOfRun * >
	Observer ()

void	slotForUpdate (const BeginOfRun *iT)

virtual	~Observer ()

Public Member Functions inherited from Observer< const BeginOfEvent * >
	Observer ()

void	slotForUpdate (const BeginOfEvent *iT)

virtual	~Observer ()

Public Member Functions inherited from Observer< const BeginOfTrack * >
	Observer ()

void	slotForUpdate (const BeginOfTrack *iT)

virtual	~Observer ()

Public Member Functions inherited from Observer< const EndOfTrack * >
	Observer ()

void	slotForUpdate (const EndOfTrack *iT)

virtual	~Observer ()

Public Member Functions inherited from Observer< const EndOfEvent * >
	Observer ()

void	slotForUpdate (const EndOfEvent *iT)

virtual	~Observer ()

Protected Member Functions
G4bool	checkHit ()

virtual void	clearHits ()

CaloG4Hit *	createNewHit ()

virtual bool	filterHit (CaloG4Hit *, double)

double	getAttenuation (G4Step *aStep, double birk1, double birk2, double birk3)

virtual uint16_t	getDepth (G4Step *)

int	getNumberOfHits ()

double	getResponseWt (G4Track *)

virtual G4bool	getStepInfo (G4Step *aStep)

virtual int	getTrackID (G4Track *)

G4bool	hitExists ()

virtual void	initRun ()

void	resetForNewPrimary (const G4ThreeVector &, double)

G4ThreeVector	setToGlobal (const G4ThreeVector &, const G4VTouchable *)

G4ThreeVector	setToLocal (const G4ThreeVector &, const G4VTouchable *)

virtual void	update (const BeginOfRun *)
	This routine will be called when the appropriate signal arrives. More...

virtual void	update (const BeginOfEvent *)
	This routine will be called when the appropriate signal arrives. More...

virtual void	update (const BeginOfTrack *trk)
	This routine will be called when the appropriate signal arrives. More...

virtual void	update (const EndOfTrack *trk)
	This routine will be called when the appropriate signal arrives. More...

virtual void	update (const ::EndOfEvent *)

void	updateHit (CaloG4Hit *)

Protected Member Functions inherited from Observer< const EndOfEvent * >
virtual void	update (const EndOfEvent *)=0
	This routine will be called when the appropriate signal arrives. More...

Protected Attributes
int	checkHits

double	correctT

bool	corrTOFBeam

CaloG4Hit *	currentHit

CaloHitID	currentID

float	edepositEM

float	edepositHAD

double	eminHit

double	eminHitD

G4int	emPDG

double	energyCut

G4ThreeVector	entranceLocal

G4ThreeVector	entrancePoint

G4int	epPDG

bool	forceSave

G4int	gammaPDG

float	incidentEnergy

double	kmaxIon

double	kmaxNeutron

double	kmaxProton

const SimTrackManager *	m_trackManager

G4ThreeVector	posGlobal

G4StepPoint *	preStepPoint

CaloHitID	previousID

int	primIDSaved

bool	runInit

bool	suppressHeavy

G4Track *	theTrack

double	tmaxHit

bool	useMap

Private Member Functions
void	cleanHitCollection ()

bool	saveHit (CaloG4Hit *)

void	storeHit (CaloG4Hit *)

void	summarize ()

Private Attributes
int	cleanIndex

int	hcID

std::map< CaloHitID, CaloG4Hit * >	hitMap

std::vector< CaloG4Hit * >	hitvec

bool	ignoreTrackID

CaloMeanResponse *	meanResponse

int	primAncestor

std::vector< CaloG4Hit * >	reusehit

std::vector< unsigned int >	selIndex

CaloSlaveSD *	slave

CaloG4HitCollection *	theHC

int	timeSlice

std::map< int, TrackWithHistory * >	tkMap

int	totalHits

Additional Inherited Members
Public Types inherited from SensitiveDetector
enum	coordinates { WorldCoordinates, LocalCoordinates }

Detailed Description

Definition at line 42 of file CaloSD.h.

Constructor & Destructor Documentation

CaloSD::CaloSD	(	G4String	aSDname,
		const DDCompactView &	cpv,
		const SensitiveDetectorCatalog &	clg,
		edm::ParameterSet const &	p,
		const SimTrackManager *	manager,
		int	tSlice = `1`,
		bool	ignoreTkID = `false`
	)

Definition at line 24 of file CaloSD.cc.

References SensitiveDetector::AssignSD(), checkHits, cleanIndex, ecaldqm::collectionName, correctT, corrTOFBeam, currentID, eminHit, eminHitD, energyCut, forceSave, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), GeV, ignoreTrackID, relval_steps::k, kmaxIon, kmaxNeutron, kmaxProton, LogDebug, SensitiveDetectorCatalog::logicalNames(), meanResponse, MeV, previousID, primAncestor, SensitiveDetector::Register(), slave, suppressHeavy, timeSlice, tmaxHit, totalHits, and useMap.

                                      : 
   SensitiveCaloDetector(name, cpv, clg, p),
   G4VGFlashSensitiveDetector(), theTrack(0), preStepPoint(0), eminHit(0), 
   eminHitD(0), m_trackManager(manager), currentHit(0), runInit(false),
   timeSlice(tSlice), ignoreTrackID(ignoreTkID), hcID(-1), theHC(0), 
   meanResponse(0) {
   //Add Hcal Sentitive Detector Names
 
   collectionName.insert(name);
 
   //Parameters
   edm::ParameterSet m_CaloSD = p.getParameter<edm::ParameterSet>("CaloSD");
   energyCut    = m_CaloSD.getParameter<double>("EminTrack")*GeV;
   tmaxHit      = m_CaloSD.getParameter<double>("TmaxHit")*ns;
   std::vector<double> eminHits = m_CaloSD.getParameter<std::vector<double> >("EminHits");
   std::vector<double> tmaxHits = m_CaloSD.getParameter<std::vector<double> >("TmaxHits");
   std::vector<std::string> hcn = m_CaloSD.getParameter<std::vector<std::string> >("HCNames");
   std::vector<int>   useResMap = m_CaloSD.getParameter<std::vector<int> >("UseResponseTables");
   std::vector<double> eminHitX = m_CaloSD.getParameter<std::vector<double> >("EminHitsDepth");
   suppressHeavy= m_CaloSD.getParameter<bool>("SuppressHeavy");
   kmaxIon      = m_CaloSD.getParameter<double>("IonThreshold")*MeV;
   kmaxProton   = m_CaloSD.getParameter<double>("ProtonThreshold")*MeV;
   kmaxNeutron  = m_CaloSD.getParameter<double>("NeutronThreshold")*MeV;
   checkHits    = m_CaloSD.getUntrackedParameter<int>("CheckHits", 25);
   useMap       = m_CaloSD.getUntrackedParameter<bool>("UseMap", true);
   int verbn    = m_CaloSD.getUntrackedParameter<int>("Verbosity", 0);
   corrTOFBeam  = m_CaloSD.getParameter<bool>("CorrectTOFBeam");
   double beamZ = m_CaloSD.getParameter<double>("BeamPosition")*cm;
   correctT     = beamZ/c_light/nanosecond;
 
   SetVerboseLevel(verbn);
   for (unsigned int k=0; k<hcn.size(); ++k) {
     if (name == (G4String)(hcn[k])) {
       if (k < eminHits.size()) eminHit = eminHits[k]*MeV;
       if (k < eminHitX.size()) eminHitD= eminHitX[k]*MeV;
       if (k < tmaxHits.size()) tmaxHit = tmaxHits[k]*ns;
       if (k < useResMap.size() && useResMap[k] > 0) meanResponse = new CaloMeanResponse(p);
       break;
     }
   }
 #ifdef DebugLog
   LogDebug("CaloSim") << "***************************************************" 
                       << "\n"
                       << "*                                                 *" 
                       << "\n"
                       << "* Constructing a CaloSD  with name " << GetName()
                       << "\n"
                       << "*                                                 *" 
                       << "\n"
                       << "***************************************************";
 #endif
   slave      = new CaloSlaveSD(name);
   currentID  = CaloHitID(timeSlice, ignoreTrackID);
   previousID = CaloHitID(timeSlice, ignoreTrackID);
   
   primAncestor = 0;
   cleanIndex = 0;
   totalHits = 0;
   forceSave = false;
 
   //
   // Now attach the right detectors (LogicalVolumes) to me
   //
   const std::vector<std::string>& lvNames = clg.logicalNames(name);
   this->Register();
   for (std::vector<std::string>::const_iterator it=lvNames.begin(); it !=lvNames.end(); ++it) {
     this->AssignSD(*it);
 #ifdef DebugLog
     LogDebug("CaloSim") << "CaloSD : Assigns SD to LV " << (*it);
 #endif
   }
 
   edm::LogInfo("CaloSim") << "CaloSD: Minimum energy of track for saving it " 
                           << energyCut/GeV  << " GeV" << "\n"
                           << "        Use of HitID Map " << useMap << "\n"
                           << "        Check last " << checkHits 
                           << " before saving the hit\n" 
                           << "        Correct TOF globally by " << correctT
                           << " ns (Flag =" << corrTOFBeam << ")\n"
                           << "        Save hits recorded before " << tmaxHit
                           << " ns and if energy is above " << eminHit/MeV
                           << " MeV (for depth 0) or " << eminHitD/MeV
                           << " MeV (for nonzero depths); Time Slice Unit " 
                           << timeSlice << " Ignore TrackID Flag " << ignoreTrackID;
 }

CaloSD::~CaloSD ( )

virtual

Definition at line 113 of file CaloSD.cc.

References meanResponse, slave, and theHC.

                 { 
   if (slave)           delete slave; 
   if (theHC)           delete theHC;
   if (meanResponse)    delete meanResponse;
 }

Member Function Documentation

G4bool CaloSD::checkHit ( )

protected

Definition at line 331 of file CaloSD.cc.

References checkHits, currentHit, currentID, newFWLiteAna::found, hitMap, j, theHC, updateHit(), and useMap.

Referenced by CastorSD::getFromLibrary(), HCalSD::getFromLibrary(), HCalSD::getFromParam(), HCalSD::getHitFibreBundle(), HCalSD::getHitPMT(), hitExists(), HCalSD::hitForFibre(), and ProcessHits().

                         {  
   //look in the HitContainer whether a hit with the same ID already exists:
   bool       found = false;
   if (useMap) {
     std::map<CaloHitID,CaloG4Hit*>::const_iterator it = hitMap.find(currentID);
     if (it != hitMap.end()) {
       currentHit = it->second;
       found      = true;
     }
   } else {
     if (checkHits <= 0) return false;
     int  minhit= (theHC->entries()>checkHits ? theHC->entries()-checkHits : 0);
     int  maxhit= theHC->entries()-1;
     
     for (int j=maxhit; j>minhit&&!found; --j) {
       if ((*theHC)[j]->getID() == currentID) {
         currentHit = (*theHC)[j];
         found      = true;
       }
     }          
   }
   
   if (found) {
     updateHit(currentHit);
     return true;
   } else {
     return false;
   }
 }

void CaloSD::cleanHitCollection ( )

private

Definition at line 688 of file CaloSD.cc.

References cleanIndex, correctT, corrTOFBeam, cond::serialization::equal(), filterHit(), CaloG4Hit::getTimeSlice(), hitMap, hitvec, i, cuy::ii, j, LogDebug, hltrates_dqm_sourceclient-live_cfg::offset, reusehit, selIndex, python.multivaluedict::sort(), swap(), theHC, cond::rpcobgas::time, totalHits, and useMap.

Referenced by EndOfEvent(), and update().

                                 {
   std::vector<CaloG4Hit*>* theCollection = theHC->GetVector();
 
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: collection before merging, size = " << theHC->entries();
 #endif
   
   selIndex.reserve(theHC->entries()-cleanIndex);
   if ( reusehit.size() == 0 ) reusehit.reserve(theHC->entries()-cleanIndex); 
 
   // if no map used, merge before hits to have the save situation as a map
   if ( !useMap ) {
     hitvec.swap(*theCollection);
     sort((hitvec.begin()+cleanIndex), hitvec.end(), CaloG4HitLess());
 #ifdef DebugLog
     LogDebug("CaloSim") << "CaloSD::cleanHitCollection: sort hits in buffer "
                         << "starting from element = " << cleanIndex;
     for (unsigned int i = 0; i<hitvec.size(); ++i) 
       LogDebug("CaloSim")<<i<<" "<<*hitvec[i];
 #endif
     unsigned int i, j;
     CaloG4HitEqual equal;
     for (i=cleanIndex; i<hitvec.size(); ++i) {
       selIndex.push_back(i-cleanIndex);
       int jump = 0;
       for (j = i+1; j <hitvec.size() && equal(hitvec[i], hitvec[j]); ++j) {
         ++jump;
         // merge j to i
         (*hitvec[i]).addEnergyDeposit(*hitvec[j]);
         (*hitvec[j]).setEM(0.);
         (*hitvec[j]).setHadr(0.);
         reusehit.push_back(hitvec[j]);
       }
       i+=jump;
     }
 #ifdef DebugLog
     LogDebug("CaloSim") << "CaloSD: cleanHitCollection merge the hits in buffer ";
     for (unsigned int i = 0; i<hitvec.size(); ++i) 
       LogDebug("CaloSim")<<i<<" "<<*hitvec[i];
 #endif
     for ( unsigned int i = cleanIndex; i < cleanIndex+selIndex.size(); ++i ) {
       hitvec[i] = hitvec[selIndex[i-cleanIndex]+cleanIndex];
     }
     hitvec.resize(cleanIndex+selIndex.size());
 #ifdef DebugLog
     LogDebug("CaloSim") << "CaloSD::cleanHitCollection: remove the merged hits in buffer,"
                         << " new size = " << hitvec.size();
     for (unsigned int i = 0; i<hitvec.size(); ++i) 
       LogDebug("CaloSim")<<i<<" "<<*hitvec[i];
 #endif
     hitvec.swap(*theCollection);
     std::vector<CaloG4Hit*>().swap(hitvec);
     selIndex.clear();
     totalHits = theHC->entries();
   }
 
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: collection after merging, size = " << theHC->entries();
 #endif
 
   int addhit = 0;
 
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: Size of reusehit after merge = " << reusehit.size();
   LogDebug("CaloSim") << "CaloSD: Starting hit selection from index = " << cleanIndex;
 #endif
   
   selIndex.reserve(theCollection->size()-cleanIndex);
   for (unsigned int i = cleanIndex; i<theCollection->size(); ++i) {   
     CaloG4Hit* aHit((*theCollection)[i]);
     
     // selection
     
     double time = aHit->getTimeSlice();
     if (corrTOFBeam) time += correctT;
     if (!filterHit(aHit,time)) {
 #ifdef DebugLog
       LogDebug("CaloSim") << "CaloSD: dropped CaloG4Hit " << " " << *aHit; 
 #endif
       
       // create the list of hits to be reused
       
       reusehit.push_back((*theCollection)[i]);
       ++addhit;
     } else {
       selIndex.push_back(i-cleanIndex);
     }
   }
 
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: Size of reusehit after selection = " << reusehit.size()  
                       << " Number of added hit = " << addhit;
 #endif
   if (useMap) {
     if ( addhit>0 ) {
       int offset = reusehit.size()-addhit;
       for (int ii = addhit-1; ii>=0; --ii) {
         CaloHitID theID = reusehit[offset+ii]->getID();
         hitMap.erase(theID);
       }
     }
   }
   for (unsigned int j = 0; j<selIndex.size(); ++j) {
     (*theCollection)[cleanIndex+j] = (*theCollection)[cleanIndex+selIndex[j]];
   }
 
   theCollection->resize(cleanIndex+selIndex.size());
   std::vector<unsigned int>().swap(selIndex);
 
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: hit collection after selection, size = "
                       << theHC->entries();
   theHC->PrintAllHits();
 #endif
     
   cleanIndex = theHC->entries();
 }

void CaloSD::clear ( void )

virtual

Definition at line 234 of file CaloSD.cc.

234 {}

void CaloSD::clearHits ( )

protectedvirtual

Implements SensitiveDetector.

Definition at line 556 of file CaloSD.cc.

References cleanIndex, hitMap, i, CaloSlaveSD::Initialize(), LogDebug, previousID, primIDSaved, CaloHitID::reset(), reusehit, slave, swap(), and useMap.

Referenced by update().

                        {  
   if (useMap) hitMap.erase (hitMap.begin(), hitMap.end());
   for (unsigned int i = 0; i<reusehit.size(); ++i) delete reusehit[i];
   std::vector<CaloG4Hit*>().swap(reusehit);
   cleanIndex  = 0;
   previousID.reset();
   primIDSaved = -99;
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: Clears hit vector for " << GetName() << " " << slave;
 #endif
   slave->Initialize();
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: Initialises slave SD for " << GetName();
 #endif
 }

CaloG4Hit * CaloSD::createNewHit ( )

protected

Definition at line 363 of file CaloSD.cc.

References currentID, CaloHitID::depth(), edepositEM, edepositHAD, energyCut, entranceLocal, entrancePoint, forceSave, GeV, incidentEnergy, LogDebug, TrackWithHistory::momentum(), NULL, posGlobal, primIDSaved, TrackInformation::putInHistory(), reusehit, TrackWithHistory::save(), CaloG4Hit::setEM(), CaloG4Hit::setEntry(), CaloG4Hit::setEntryLocal(), CaloG4Hit::setHadr(), CaloG4Hit::setID(), CaloG4Hit::setIncidentEnergy(), CaloG4Hit::setPosition(), mathSSE::sqrt(), storeHit(), TrackInformation::storeTrack(), theTrack, CaloHitID::timeSliceID(), tkMap, totalHits, CaloHitID::trackID(), CaloHitID::unitID(), updateHit(), and useMap.

Referenced by ZdcSD::getFromLibrary(), CastorSD::getFromLibrary(), HCalSD::getFromLibrary(), HCalSD::getFromParam(), HCalSD::getHitFibreBundle(), HCalSD::getHitPMT(), HCalSD::hitForFibre(), ZdcSD::ProcessHits(), DreamSD::ProcessHits(), HCalSD::ProcessHits(), and ProcessHits().

                                 {
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD::CreateNewHit for"
                       << " Unit " << currentID.unitID() 
                       << " " << currentID.depth()
                       << " Edeposit = " << edepositEM << " " << edepositHAD;
   LogDebug("CaloSim") << " primary "    << currentID.trackID()
                       << " time slice " << currentID.timeSliceID()
                       << " For Track  " << theTrack->GetTrackID()
                       << " which is a " <<theTrack->GetDefinition()->GetParticleName()
                       << " of energy "  << theTrack->GetKineticEnergy()/GeV
                       << " " << theTrack->GetMomentum().mag()/GeV
                       << " daughter of part. " << theTrack->GetParentID()
                       << " and created by " ;
   
   if (theTrack->GetCreatorProcess()!=NULL)
     LogDebug("CaloSim") << theTrack->GetCreatorProcess()->GetProcessName() ;
   else 
     LogDebug("CaloSim") << "NO process";
 #endif  
   
   CaloG4Hit* aHit;
   if (reusehit.size() > 0) {
     aHit = reusehit[0];
     aHit->setEM(0.);
     aHit->setHadr(0.);
     reusehit.erase(reusehit.begin());
   } else {
     aHit = new CaloG4Hit;
   }
   
   aHit->setID(currentID);
   aHit->setEntry(entrancePoint.x(),entrancePoint.y(),entrancePoint.z());
   aHit->setEntryLocal(entranceLocal.x(),entranceLocal.y(),entranceLocal.z());
   aHit->setPosition(posGlobal.x(),posGlobal.y(),posGlobal.z());
   aHit->setIncidentEnergy(incidentEnergy);
   updateHit(aHit);
   
   storeHit(aHit);
   double etrack = 0;
   if (currentID.trackID() == primIDSaved) { // The track is saved; nothing to be done
   } else if (currentID.trackID() == theTrack->GetTrackID()) {
     etrack= theTrack->GetKineticEnergy();
     //edm::LogInfo("CaloSim") << "CaloSD: set save the track " << currentID.trackID()
     //      << " etrack " << etrack << " eCut " << energyCut << " flag " << forceSave;
     if (etrack >= energyCut || forceSave) {
       TrackInformation* trkInfo = (TrackInformation *)(theTrack->GetUserInformation());
       trkInfo->storeTrack(true);
       trkInfo->putInHistory();
       //      trkInfo->setAncestor();
 #ifdef DebugLog
       LogDebug("CaloSim") << "CaloSD: set save the track " << currentID.trackID()
                           << " with Hit";
 #endif
     }
   } else {
     TrackWithHistory * trkh = tkMap[currentID.trackID()];
 #ifdef DebugLog
     LogDebug("CaloSim") << "CaloSD : TrackwithHistory pointer for " 
                         << currentID.trackID() << " is " << trkh;
 #endif
     if (trkh != NULL) {
       etrack = sqrt(trkh->momentum().Mag2());
       if (etrack >= energyCut) {
         trkh->save();
 #ifdef DebugLog
         LogDebug("CaloSim") << "CaloSD: set save the track " 
                             << currentID.trackID() << " with Hit";
 #endif
       }
     }
   }
   primIDSaved = currentID.trackID();
   if (useMap) totalHits++;
   return aHit;
 }  

void CaloSD::DrawAll ( )

virtual

Definition at line 236 of file CaloSD.cc.

236 {}

void CaloSD::EndOfEvent ( G4HCofThisEvent * eventHC )

virtual

Reimplemented from SensitiveDetector.

Definition at line 222 of file CaloSD.cc.

References cleanHitCollection(), and theHC.

                                          {
   // clean the hits for the last tracks
   
   cleanHitCollection();
   
 #ifdef DebugLog
   edm::LogInfo("CaloSim") << "CaloSD: EndofEvent entered with " << theHC->entries()
                           << " entries";
 #endif
   //  TimeMe("CaloSD:sortAndMergeHits",false);
 }

void CaloSD::fillHits	(	edm::PCaloHitContainer &	c,
		std::string	n
	)

virtual

Implements SensitiveCaloDetector.

Definition at line 245 of file CaloSD.cc.

References CaloSlaveSD::Clean(), CaloSlaveSD::hits(), gen::n, CaloSlaveSD::name(), and slave.

                                                           {
   if (slave->name() == n) c=slave->hits();
   slave->Clean();
 }

bool CaloSD::filterHit	(	CaloG4Hit *	hit,
		double	time
	)

protectedvirtual

Reimplemented in HCalSD.

Definition at line 597 of file CaloSD.cc.

References eminHit, eminHitD, CaloG4Hit::getDepth(), CaloG4Hit::getEnergyDeposit(), LogDebug, and tmaxHit.

Referenced by cleanHitCollection().

                                                   {
   double emin(eminHit);
   if (hit->getDepth() > 0) emin = eminHitD;
 #ifdef DebugLog
   LogDebug("CaloSim") << "Depth " << hit->getDepth() << " Emin = " << emin << " ("
                       << eminHit << ", " << eminHitD << ")";
 #endif   
   return ((time <= tmaxHit) && (hit->getEnergyDeposit() > emin));
 }

double CaloSD::getAttenuation	(	G4Step *	aStep,
		double	birk1,
		double	birk2,
		double	birk3
	)

protected

Definition at line 465 of file CaloSD.cc.

References funct::abs(), EnergyCorrector::c, RecoTauCleanerPlugins::charge, LogDebug, and histoStyle::weight.

Referenced by EcalTBH4BeamSD::getEnergyDeposit(), HcalTB06BeamSD::getEnergyDeposit(), ECalSD::getEnergyDeposit(), HcalTB02SD::getEnergyDeposit(), HCalSD::getEnergyDeposit(), and DreamSD::getStepInfo().

                                                                                      {
   double weight = 1.;
   double charge = aStep->GetPreStepPoint()->GetCharge();
 
   if (charge != 0. && aStep->GetStepLength() > 0) {
     G4Material* mat = aStep->GetPreStepPoint()->GetMaterial();
     double density = mat->GetDensity();
     double dedx    = aStep->GetTotalEnergyDeposit()/aStep->GetStepLength();
     double rkb     = birk1/density;
     double c       = birk2*rkb*rkb;
     if (std::abs(charge) >= 2.) rkb /= birk3; // based on alpha particle data
     weight = 1./(1.+rkb*dedx+c*dedx*dedx);
 #ifdef DebugLog
     LogDebug("CaloSim") << "CaloSD::getAttenuation in " << mat->GetName() 
                         << " Charge " << charge << " dE/dx " << dedx 
                         << " Birk Const " << rkb << ", " << c << " Weight = " 
                         << weight << " dE " << aStep->GetTotalEnergyDeposit();
 #endif
   }
   return weight;
 }

uint16_t CaloSD::getDepth ( G4Step * )

protectedvirtual

Reimplemented in ECalSD.

Definition at line 595 of file CaloSD.cc.

Referenced by getStepInfo(), and ProcessHits().

595 { return 0; }

double CaloSD::getEnergyDeposit ( G4Step * step )

virtual

Reimplemented in HCalSD, CastorSD, HcalTB02SD, ECalSD, EcalTBH4BeamSD, and HcalTB06BeamSD.

Definition at line 203 of file CaloSD.cc.

Referenced by getStepInfo().

                                              {
   return aStep->GetTotalEnergyDeposit();
 }

int CaloSD::getNumberOfHits ( )

protected

Definition at line 361 of file CaloSD.cc.

References theHC.

Referenced by HCalSD::ProcessHits().

361 { return theHC->entries(); }

CaloSD::theHC

CaloG4HitCollection * theHC

Definition: CaloSD.h:145

double CaloSD::getResponseWt ( G4Track * aTrack )

protected

Definition at line 607 of file CaloSD.cc.

References TrackInformation::genParticleP(), TrackInformation::genParticlePID(), CaloMeanResponse::getWeight(), and meanResponse.

Referenced by ECalSD::getEnergyDeposit(), and HCalSD::getEnergyDeposit().

                                             {
   if (meanResponse) {
     TrackInformation * trkInfo = (TrackInformation *)(aTrack->GetUserInformation());
     return meanResponse->getWeight(trkInfo->genParticlePID(), trkInfo->genParticleP());
   } else {
     return 1;
   }
 }

bool CaloSD::getStepInfo ( G4Step * aStep )

protectedvirtual

Reimplemented in DreamSD.

Definition at line 250 of file CaloSD.cc.

References currentID, TauDecayModes::dec, HLT_25ns14e33_v1_cff::depth, edepositEM, edepositHAD, emPDG, epPDG, gammaPDG, getDepth(), getEnergyDeposit(), getTrackID(), LogDebug, preStepPoint, setDetUnitId(), CaloHitID::setID(), theTrack, cond::rpcobgas::time, and CaloHitID::unitID().

Referenced by ZdcSD::ProcessHits(), HCalSD::ProcessHits(), and ProcessHits().

                                       {  
 
   preStepPoint = aStep->GetPreStepPoint(); 
   theTrack     = aStep->GetTrack();   
   
   double       time  = (aStep->GetPostStepPoint()->GetGlobalTime())/nanosecond;
   unsigned int unitID= setDetUnitId(aStep);
   uint16_t     depth = getDepth(aStep);
   int          primaryID = getTrackID(theTrack);
   
   bool flag = (unitID > 0);
   if (flag) {
     currentID.setID(unitID, time, primaryID, depth);
 #ifdef DebugLog
     G4TouchableHistory* touch =(G4TouchableHistory*)(theTrack->GetTouchable());
     LogDebug("CaloSim") << "CaloSD:: GetStepInfo for"
                         << " PV "     << touch->GetVolume(0)->GetName()
                         << " PVid = " << touch->GetReplicaNumber(0)
                         << " MVid = " << touch->GetReplicaNumber(1)
                         << " Unit   " << currentID.unitID() 
                         << " Edeposit = " << edepositEM << " " << edepositHAD;
   } else {
     G4TouchableHistory* touch =(G4TouchableHistory*)(theTrack->GetTouchable());
     LogDebug("CaloSim") << "CaloSD:: GetStepInfo for"
                         << " PV "     << touch->GetVolume(0)->GetName()
                         << " PVid = " << touch->GetReplicaNumber(0)
                         << " MVid = " << touch->GetReplicaNumber(1)
                         << " Unit   " << std::hex << unitID << std::dec 
                         << " Edeposit = " << edepositEM << " " << edepositHAD;
 #endif
   }
   
   G4int particleCode = theTrack->GetDefinition()->GetPDGEncoding();
   if (particleCode == emPDG ||
       particleCode == epPDG ||
       particleCode == gammaPDG ) {
     edepositEM  = getEnergyDeposit(aStep);
     edepositHAD = 0.;
   } else {
     edepositEM  = 0.;
     edepositHAD = getEnergyDeposit(aStep);
   }
 
   return flag;
 }

int CaloSD::getTrackID ( G4Track * aTrack )

protectedvirtual

Reimplemented in ECalSD.

Definition at line 574 of file CaloSD.cc.

References forceSave, TrackInformation::getIDonCaloSurface(), LogDebug, and preStepPoint.

Referenced by getStepInfo(), ECalSD::getTrackID(), and ProcessHits().

                                       {
   int primaryID = 0;
   forceSave = false;
   TrackInformation* trkInfo=(TrackInformation *)(aTrack->GetUserInformation());
   if (trkInfo) {
     primaryID = trkInfo->getIDonCaloSurface(); 
 #ifdef DebugLog
     LogDebug("CaloSim") << "CaloSD: hit update from track Id on Calo Surface " 
                         << trkInfo->getIDonCaloSurface();
 #endif   
   } else {
     primaryID = aTrack->GetTrackID();
 #ifdef DebugLog
     edm::LogWarning("CaloSim") << "CaloSD: Problem with primaryID **** set by "
                                << "force to TkID **** " << primaryID << " in "
                                << preStepPoint->GetTouchable()->GetVolume(0)->GetName();
 #endif
   }
   return primaryID;
 }

G4bool CaloSD::hitExists ( )

protected

Definition at line 310 of file CaloSD.cc.

References checkHit(), currentHit, currentID, posGlobal, preStepPoint, previousID, resetForNewPrimary(), CaloHitID::trackID(), and updateHit().

Referenced by ZdcSD::ProcessHits(), DreamSD::ProcessHits(), HCalSD::ProcessHits(), and ProcessHits().

                          {
 #ifdef DebugLog
   if (currentID.trackID()<1)
     edm::LogWarning("CaloSim") << "***** CaloSD error: primaryID = " 
                                << currentID.trackID()
                                << " maybe detector name changed";
 #endif  
   // Update if in the same detector, time-slice and for same track   
   if (currentID == previousID) {
     updateHit(currentHit);
     return true;
   }
   
   // Reset entry point for new primary
   posGlobal = preStepPoint->GetPosition();
   if (currentID.trackID() != previousID.trackID()) 
     resetForNewPrimary(preStepPoint->GetPosition(), preStepPoint->GetKineticEnergy());
   
   return checkHit();
 }

void CaloSD::Initialize ( G4HCofThisEvent * HCE )

virtual

Reimplemented from SensitiveDetector.

Definition at line 207 of file CaloSD.cc.

References ecaldqm::collectionName, hcID, theHC, and totalHits.

                                              { 
   totalHits = 0;
   
 #ifdef DebugLog
   edm::LogInfo("CaloSim") << "CaloSD : Initialize called for " << GetName(); 
 #endif
   
   //This initialization is performed at the beginning of an event
   //------------------------------------------------------------
   theHC = new CaloG4HitCollection(GetName(), collectionName[0]);
   
   if (hcID<0) hcID = G4SDManager::GetSDMpointer()->GetCollectionID(collectionName[0]);
   HCE->AddHitsCollection(hcID, theHC);
 }

void CaloSD::initRun ( void )

protectedvirtual

Reimplemented in CastorSD, HCalSD, DreamSD, and ZdcSD.

Definition at line 572 of file CaloSD.cc.

Referenced by update().

572 {}

void CaloSD::PrintAll ( )

virtual

Definition at line 238 of file CaloSD.cc.

References theHC.

                       {
 #ifdef DebugLog
   edm::LogInfo("CaloSim") << "CaloSD: Collection " << theHC->GetName();
 #endif
   theHC->PrintAllHits();
 } 

bool CaloSD::ProcessHits	(	G4Step *	step,
		G4TouchableHistory *	tHistory
	)

virtual

Implements SensitiveDetector.

Reimplemented in HCalSD, DreamSD, and ZdcSD.

Definition at line 119 of file CaloSD.cc.

References createNewHit(), currentHit, edepositEM, edepositHAD, getStepInfo(), hitExists(), SensitiveDetector::NaNTrap(), and NULL.

                                                               {
   
   NaNTrap( aStep ) ;
   
   if (aStep == NULL) {
     return true;
   } else {
     if (getStepInfo(aStep)) {
       if (hitExists() == false && edepositEM+edepositHAD>0.) 
         currentHit = createNewHit();
     }
   }
   return true;
 } 

bool CaloSD::ProcessHits	(	G4GFlashSpot *	aSpot,
		G4TouchableHistory *
	)

virtual

Definition at line 134 of file CaloSD.cc.

References checkHit(), createNewHit(), currentHit, currentID, TauDecayModes::dec, HLT_25ns14e33_v1_cff::depth, edepositEM, edepositHAD, emPDG, entranceLocal, entrancePoint, epPDG, gammaPDG, getDepth(), getTrackID(), GeV, incidentEnergy, LogDebug, NULL, posGlobal, preStepPoint, previousID, setDetUnitId(), CaloHitID::setID(), theTrack, cond::rpcobgas::time, CaloHitID::trackID(), CaloHitID::unitID(), and updateHit().

                                                                  { 
 
   if (aSpot != NULL) {   
     theTrack = const_cast<G4Track *>(aSpot->GetOriginatorTrack()->GetPrimaryTrack());
     G4int particleCode = theTrack->GetDefinition()->GetPDGEncoding();
     
     if (particleCode == emPDG ||
         particleCode == epPDG ||
         particleCode == gammaPDG ) {
       edepositEM  = aSpot->GetEnergySpot()->GetEnergy();
       edepositHAD = 0.;
     } else {
       edepositEM  = 0.;
       edepositHAD = 0.;
     }
  
     if (edepositEM>0.) {
       G4Step *      fFakeStep          = new G4Step();
       preStepPoint                     = fFakeStep->GetPreStepPoint();
       G4StepPoint * fFakePostStepPoint = fFakeStep->GetPostStepPoint();
       preStepPoint->SetPosition(aSpot->GetPosition());
       fFakePostStepPoint->SetPosition(aSpot->GetPosition());
       
       G4TouchableHandle fTouchableHandle   = aSpot->GetTouchableHandle();
       preStepPoint->SetTouchableHandle(fTouchableHandle);
       fFakeStep->SetTotalEnergyDeposit(aSpot->GetEnergySpot()->GetEnergy());
       
       double       time   = 0;
       unsigned int unitID = setDetUnitId(fFakeStep);
       int          primaryID = getTrackID(theTrack);
       uint16_t     depth = getDepth(fFakeStep);
 
       if (unitID > 0) {
         currentID.setID(unitID, time, primaryID, depth);
 #ifdef DebugLog
         LogDebug("CaloSim") << "CaloSD:: GetSpotInfo for"
                             << " Unit 0x" << std::hex << currentID.unitID() 
                             << std::dec << " Edeposit = " << edepositEM << " " 
                             << edepositHAD;
 #endif
         // Update if in the same detector, time-slice and for same track   
         if (currentID == previousID) {
       updateHit(currentHit);
         } else {
           posGlobal = aSpot->GetEnergySpot()->GetPosition();
           // Reset entry point for new primary
           if (currentID.trackID() != previousID.trackID()) {
             entrancePoint  = aSpot->GetPosition();
             entranceLocal  = aSpot->GetTouchableHandle()->GetHistory()->
                                       GetTopTransform().TransformPoint(entrancePoint);
             incidentEnergy = theTrack->GetKineticEnergy();
 #ifdef DebugLog
             LogDebug("CaloSim") << "CaloSD: Incident energy " 
                                 << incidentEnergy/GeV << " GeV and" 
                                 << " entrance point " << entrancePoint 
                                 << " (Global) " << entranceLocal << " (Local)";
 #endif
           }
 
           if (checkHit() == false) currentHit = createNewHit();
         }
       }
       delete  fFakeStep;
     }
     return true;
   } 
   return false;
 }                                   

void CaloSD::resetForNewPrimary	(	const G4ThreeVector &	point,
		double	energy
	)

protected

Definition at line 454 of file CaloSD.cc.

References relval_parameters_module::energy, entranceLocal, entrancePoint, GeV, incidentEnergy, LogDebug, point, preStepPoint, and setToLocal().

Referenced by ZdcSD::getFromLibrary(), CastorSD::getFromLibrary(), HCalSD::getFromLibrary(), HCalSD::getHitFibreBundle(), HCalSD::getHitPMT(), hitExists(), ZdcSD::setTrackID(), CastorSD::setTrackID(), and HCalSD::setTrackID().

                                                                          { 
   entrancePoint  = point;
   entranceLocal  = setToLocal(entrancePoint, preStepPoint->GetTouchable());
   incidentEnergy = energy;
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: Incident energy " << incidentEnergy/GeV 
                       << " GeV and" << " entrance point " << entrancePoint 
                       << " (Global) " << entranceLocal << " (Local)";
 #endif
 }

bool CaloSD::saveHit ( CaloG4Hit * aHit )

private

Definition at line 627 of file CaloSD.cc.

References correctT, corrTOFBeam, TauDecayModes::dec, CaloG4Hit::getDepth(), CaloG4Hit::getEM(), CaloG4Hit::getHadr(), CaloG4Hit::getTimeSlice(), CaloG4Hit::getTrackID(), CaloG4Hit::getUnitID(), GeV, SimTrackManager::giveMotherNeeded(), LogDebug, m_trackManager, convertSQLiteXML::ok, CaloSlaveSD::processHits(), slave, cond::rpcobgas::time, and SimTrackManager::trackExists().

Referenced by update().

                                     {  
   int tkID;
   bool ok   = true;
   if (m_trackManager) {
     tkID = m_trackManager->giveMotherNeeded(aHit->getTrackID());
     if (tkID == 0) {
       if (m_trackManager->trackExists(aHit->getTrackID())) tkID = (aHit->getTrackID());
       else {
     ok = false;
       }
     }
   } else {
     tkID = aHit->getTrackID();
     ok = false;
   }
   //  edm::LogInfo("CaloSim") << "CalosD: Track ID " << aHit->getTrackID() << " changed to " << tkID << " by SimTrackManager" << " Status " << ok;
 #ifdef DebugLog
   LogDebug("CaloSim") << "CalosD: Track ID " << aHit->getTrackID() 
                       << " changed to " << tkID << " by SimTrackManager"
                       << " Status " << ok;
 #endif
   double time = aHit->getTimeSlice();
   if (corrTOFBeam) time += correctT;
   slave->processHits(aHit->getUnitID(), aHit->getEM()/GeV, 
                      aHit->getHadr()/GeV, time, tkID, aHit->getDepth());
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: Store Hit at " << std::hex 
                       << aHit->getUnitID() << std::dec << " " 
                       << aHit->getDepth() << " due to " << tkID 
                       << " in time " << time << " of energy " 
                       << aHit->getEM()/GeV << " GeV (EM) and " 
                       << aHit->getHadr()/GeV << " GeV (Hadr)";
 #endif
   return ok;
 }

virtual uint32_t CaloSD::setDetUnitId ( G4Step * step )

pure virtual

Implements SensitiveDetector.

Implemented in HCalSD, CastorSD, HcalTB02SD, ECalSD, EcalTBH4BeamSD, HcalTB06BeamSD, DreamSD, and ZdcSD.

Referenced by getStepInfo(), and ProcessHits().

G4ThreeVector CaloSD::setToGlobal	(	const G4ThreeVector &	local,
		const G4VTouchable *	touch
	)

protected

Definition at line 303 of file CaloSD.cc.

                                                                                        {
 
   G4ThreeVector globalPoint = touch->GetHistory()->GetTopTransform().Inverse().TransformPoint(local);
   
   return globalPoint;  
 }

G4ThreeVector CaloSD::setToLocal	(	const G4ThreeVector &	global,
		const G4VTouchable *	touch
	)

protected

Definition at line 296 of file CaloSD.cc.

Referenced by DreamSD::curve_LY(), HcalTB02SD::curve_LY(), ECalSD::curve_LY(), ECalSD::getRadiationLength(), resetForNewPrimary(), and HcalTB06BeamSD::setDetUnitId().

                                                                                        {
 
   G4ThreeVector localPoint = touch->GetHistory()->GetTopTransform().TransformPoint(global);
   
   return localPoint;  
 }

void CaloSD::storeHit ( CaloG4Hit * hit )

private

Definition at line 616 of file CaloSD.cc.

References hitMap, previousID, theHC, CaloHitID::trackID(), and useMap.

Referenced by createNewHit().

                                     {
   if (previousID.trackID()<0) return;
   if (hit == 0) {
     edm::LogWarning("CaloSim") << "CaloSD: hit to be stored is NULL !!";
     return;
   }
   
   theHC->insert(hit);
   if (useMap) hitMap.insert(std::pair<CaloHitID,CaloG4Hit*>(previousID,hit));
 }

void CaloSD::summarize ( )

private

Definition at line 663 of file CaloSD.cc.

Referenced by update().

663 {}

void CaloSD::update ( const BeginOfRun * )

protectedvirtual

This routine will be called when the appropriate signal arrives.

Implements Observer< const BeginOfRun * >.

Definition at line 487 of file CaloSD.cc.

References emPDG, epPDG, gammaPDG, initRun(), LogDebug, and runInit.

                                       {
   G4ParticleTable * theParticleTable = G4ParticleTable::GetParticleTable();
   G4String particleName;
   emPDG = theParticleTable->FindParticle(particleName="e-")->GetPDGEncoding();
   epPDG = theParticleTable->FindParticle(particleName="e+")->GetPDGEncoding();
   gammaPDG = theParticleTable->FindParticle(particleName="gamma")->GetPDGEncoding();
 #ifdef DebugLog
   LogDebug("CaloSim") << "CaloSD: Particle code for e- = " << emPDG
                       << " for e+ = " << epPDG << " for gamma = " << gammaPDG;
 #endif
   initRun();
   runInit = true;
 } 

void CaloSD::update ( const BeginOfEvent * )

protectedvirtual

This routine will be called when the appropriate signal arrives.

Implements Observer< const BeginOfEvent * >.

Definition at line 501 of file CaloSD.cc.

References clearHits(), and LogDebug.

                                         {
 #ifdef DebugLog
   LogDebug("CaloSim")  << "CaloSD: Dispatched BeginOfEvent for " << GetName() 
                        << " !" ;
 #endif
   clearHits();
 }

void CaloSD::update ( const BeginOfTrack * )

protectedvirtual

This routine will be called when the appropriate signal arrives.

Implements Observer< const BeginOfTrack * >.

Definition at line 665 of file CaloSD.cc.

References cleanHitCollection(), TrackInformation::isPrimary(), LogDebug, primAncestor, and theHC.

                                             {
   int primary = -1;
   TrackInformation * trkInfo = (TrackInformation *)((*trk)()->GetUserInformation());
   if ( trkInfo->isPrimary() ) primary = (*trk)()->GetTrackID();
   
 #ifdef DebugLog
   LogDebug("CaloSim") << "New track: isPrimary " << trkInfo->isPrimary() 
                       << " primary ID = " << primary 
                       << " primary ancestor ID " << primAncestor;
 #endif
   
   // update the information if a different primary track ID 
   
   if (primary > 0 && primary != primAncestor) {
     primAncestor = primary;
     
     // clean the hits information
     
     if (theHC->entries()>0) cleanHitCollection();
     
   }
 }

void CaloSD::update ( const EndOfTrack * )

protectedvirtual

This routine will be called when the appropriate signal arrives.

Implements Observer< const EndOfTrack * >.

Definition at line 509 of file CaloSD.cc.

References TrackInformation::getIDonCaloSurface(), LogDebug, m_trackManager, NULL, tkMap, SimTrackManager::trackContainer(), and TrackWithHistory::trackID().

                                           {
   int id = (*trk)()->GetTrackID();
   TrackInformation *trkI =(TrackInformation *)((*trk)()->GetUserInformation());
   int lastTrackID = -1;
   if (trkI) lastTrackID = trkI->getIDonCaloSurface();
   if (id == lastTrackID) {
     const TrackContainer * trksForThisEvent = m_trackManager->trackContainer();
     if (trksForThisEvent != NULL) {
       int it = (int)(trksForThisEvent->size()) - 1;
       if (it >= 0) {
         TrackWithHistory * trkH = (*trksForThisEvent)[it];
         if (trkH->trackID() == (unsigned int)(id)) tkMap[id] = trkH;
 #ifdef DebugLog
         LogDebug("CaloSim") << "CaloSD: get track " << it << " from "
                             << "Container of size " << trksForThisEvent->size()
                             << " with ID " << trkH->trackID();
       } else {
         LogDebug("CaloSim") << "CaloSD: get track " << it << " from "
                             << "Container of size " << trksForThisEvent->size()
                             << " with no ID";
 #endif
       }
     }
   }
 }

void CaloSD::update ( const ::EndOfEvent * )

protectedvirtual

Definition at line 535 of file CaloSD.cc.

References prof2calltree::count, i, convertSQLiteXML::ok, CaloSlaveSD::ReserveMemory(), saveHit(), slave, summarize(), theHC, tkMap, and totalHits.

                                        {
   int count = 0, wrong = 0;
   bool ok;
   
   slave->ReserveMemory(theHC->entries());
 
   for (int i=0; i<theHC->entries(); ++i) {
     ok = saveHit((*theHC)[i]);
     ++count;
     if (!ok)  ++wrong;
   }
   
   edm::LogInfo("CaloSim") << "CaloSD: " << GetName() << " store " << count
                           << " hits recorded with " << wrong 
                           << " track IDs not given properly and "
                           << totalHits-count << " hits not passing cuts";
   summarize();
 
   tkMap.erase (tkMap.begin(), tkMap.end());
 }

void CaloSD::updateHit ( CaloG4Hit * aHit )

protected

Definition at line 440 of file CaloSD.cc.

References CaloG4Hit::addEnergyDeposit(), currentID, edepositEM, edepositHAD, LogDebug, MeV, and previousID.

Referenced by checkHit(), createNewHit(), ZdcSD::getFromLibrary(), CastorSD::getFromLibrary(), HCalSD::getFromLibrary(), HCalSD::getFromParam(), HCalSD::getHitFibreBundle(), HCalSD::getHitPMT(), hitExists(), HCalSD::hitForFibre(), and ProcessHits().

                                       {
   if (edepositEM+edepositHAD != 0) {
     aHit->addEnergyDeposit(edepositEM,edepositHAD);
 #ifdef DebugLog
     LogDebug("CaloSim") << "CaloSD: Add energy deposit in " << currentID 
                         << " em " << edepositEM/MeV << " hadronic " 
                         << edepositHAD/MeV << " MeV"; 
 #endif
   }
 
   // buffer for next steps:
   previousID = currentID;
 }