|
|
|
#include <SimG4CMS/Forward/interface/CastorSD.h>
Public Member Functions | |
| CastorSD (G4String, const DDCompactView &, SensitiveDetectorCatalog &clg, edm::ParameterSet const &, const SimTrackManager *) | |
| virtual double | getEnergyDeposit (G4Step *) |
| virtual uint32_t | setDetUnitId (G4Step *step) |
| void | setNumberingScheme (CastorNumberingScheme *scheme) |
| virtual | ~CastorSD () |
Protected Member Functions | |
| virtual void | initRun () |
Private Member Functions | |
| void | getFromLibrary (G4Step *) |
| uint32_t | rotateUnitID (uint32_t, G4Track *, CastorShowerEvent) |
| int | setTrackID (G4Step *) |
Private Attributes | |
| double | energyThresholdSL |
| G4LogicalVolume * | lvC3EF |
| G4LogicalVolume * | lvC3HF |
| G4LogicalVolume * | lvC4EF |
| G4LogicalVolume * | lvC4HF |
| G4LogicalVolume * | lvCAST |
| double | non_compensation_factor |
| CastorNumberingScheme * | numberingScheme |
| CastorShowerLibrary * | showerLibrary |
| bool | useShowerLibrary |
Description: Stores hits of Castor in appropriate container
Usage: Used in sensitive detector builder
Definition at line 30 of file CastorSD.h.
| CastorSD::CastorSD | ( | G4String | name, |
| const DDCompactView & | cpv, | ||
| SensitiveDetectorCatalog & | clg, | ||
| edm::ParameterSet const & | p, | ||
| const SimTrackManager * | manager | ||
| ) |
Definition at line 30 of file CastorSD.cc.
References energyThresholdSL, edm::ParameterSet::getParameter(), lvC3EF, lvC3HF, lvC4EF, lvC4HF, lvCAST, non_compensation_factor, setNumberingScheme(), showerLibrary, and useShowerLibrary.
: CaloSD(name, cpv, clg, p, manager), numberingScheme(0), lvC3EF(0), lvC3HF(0), lvC4EF(0), lvC4HF(0), lvCAST(0) { edm::ParameterSet m_CastorSD = p.getParameter<edm::ParameterSet>("CastorSD"); useShowerLibrary = m_CastorSD.getParameter<bool>("useShowerLibrary"); energyThresholdSL = m_CastorSD.getParameter<double>("minEnergyInGeVforUsingSLibrary"); energyThresholdSL = energyThresholdSL*GeV; // Convert GeV => MeV non_compensation_factor = m_CastorSD.getParameter<double>("nonCompensationFactor"); if (useShowerLibrary) showerLibrary = new CastorShowerLibrary(name, p); setNumberingScheme(new CastorNumberingScheme()); edm::LogInfo("ForwardSim") << "***************************************************\n" << "* *\n" << "* Constructing a CastorSD with name " << GetName() << "\n" << "* *\n" << "***************************************************"; const G4LogicalVolumeStore * lvs = G4LogicalVolumeStore::GetInstance(); std::vector<G4LogicalVolume*>::const_iterator lvcite; for (lvcite = lvs->begin(); lvcite != lvs->end(); lvcite++) { if (strcmp(((*lvcite)->GetName()).c_str(),"C3EF") == 0) lvC3EF = (*lvcite); if (strcmp(((*lvcite)->GetName()).c_str(),"C3HF") == 0) lvC3HF = (*lvcite); if (strcmp(((*lvcite)->GetName()).c_str(),"C4EF") == 0) lvC4EF = (*lvcite); if (strcmp(((*lvcite)->GetName()).c_str(),"C4HF") == 0) lvC4HF = (*lvcite); if (strcmp(((*lvcite)->GetName()).c_str(),"CAST") == 0) lvCAST = (*lvcite); if (lvC3EF != 0 && lvC3HF != 0 && lvC4EF != 0 && lvC4HF != 0 && lvCAST != 0) break; } edm::LogInfo("ForwardSim") << "CastorSD:: LogicalVolume pointers\n" << lvC3EF << " for C3EF; " << lvC3HF << " for C3HF; " << lvC4EF << " for C4EF; " << lvC4HF << " for C4HF; " << lvCAST << " for CAST. " << std::endl; // if(useShowerLibrary) edm::LogInfo("ForwardSim") << "\n Using Castor Shower Library \n"; }
| CastorSD::~CastorSD | ( | ) | [virtual] |
Definition at line 77 of file CastorSD.cc.
References showerLibrary, and useShowerLibrary.
{
if (useShowerLibrary) delete showerLibrary;
}
| double CastorSD::getEnergyDeposit | ( | G4Step * | aStep | ) | [virtual] |
Reimplemented from CaloSD.
Definition at line 94 of file CastorSD.cc.
References a, abs, beta, DeDxDiscriminatorTools::charge(), dot(), energyThresholdSL, eta(), getFromLibrary(), funct::log(), LogDebug, lvC3EF, lvC3HF, lvC4EF, lvC4HF, lvCAST, M_PI, max(), min, SensitiveDetector::name, NULL, phi, pi, CaloSD::preStepPoint, dttmaxenums::R, alignCSCRings::r, mathSSE::sqrt(), funct::tan(), theta(), CaloSD::theTrack, and useShowerLibrary.
{
float NCherPhot = 0.;
if (aStep == NULL) {
return 0;
} else {
// preStepPoint information *********************************************
G4StepPoint* preStepPoint = aStep->GetPreStepPoint();
G4VPhysicalVolume* currentPV = preStepPoint->GetPhysicalVolume();
G4LogicalVolume* currentLV = currentPV->GetLogicalVolume();
G4String name = currentPV->GetName();
std::string nameVolume;
nameVolume.assign(name,0,4);
#ifdef debugLog
G4SteppingControl stepControlFlag = aStep->GetControlFlag();
if (aStep->IsFirstStepInVolume())
LogDebug("ForwardSim") << "CastorSD::getEnergyDeposit:"
<< "\n IsFirstStepInVolume " ;
#endif
// Get theTrack
G4Track* theTrack = aStep->GetTrack();
#ifdef debugLog
if (useShowerLibrary && currentLV==lvCAST) {
LogDebug("ForwardSim") << "CastorSD::getEnergyDeposit:"
<< "\n TrackID , ParentID , ParticleName ,"
<< " eta , phi , z , time ,"
<< " K , E , Mom "
<< "\n TRACKINFO: "
<< theTrack->GetTrackID()
<< " , "
<< theTrack->GetParentID()
<< " , "
<< theTrack->GetDefinition()->GetParticleName()
<< " , "
<< theTrack->GetPosition().eta()
<< " , "
<< theTrack->GetPosition().phi()
<< " , "
<< theTrack->GetPosition().z()
<< " , "
<< theTrack->GetGlobalTime()
<< " , "
<< theTrack->GetKineticEnergy()
<< " , "
<< theTrack->GetTotalEnergy()
<< " , "
<< theTrack->GetMomentum().mag() ;
if(theTrack->GetTrackID() != 1)
LogDebug("ForwardSim") << "CastorSD::getEnergyDeposit:"
<< "\n CurrentStepNumber , TrackID , Particle , VertexPosition ,"
<< " LogicalVolumeAtVertex , CreatorProcess"
<< "\n TRACKINFO2: "
<< theTrack->GetCurrentStepNumber()
<< " , "
<< theTrack->GetTrackID()
<< " , "
<< theTrack->GetDefinition()->GetParticleName()
<< " , "
<< theTrack->GetVertexPosition()
<< " , "
<< theTrack->GetLogicalVolumeAtVertex()->GetName()
<< " , "
<< theTrack->GetCreatorProcess()->GetProcessName() ;
} // end of if(useShowerLibrary)
#endif
// if particle moves from interaction point or "backwards (halo)
bool backward = false;
G4ThreeVector hitPoint = preStepPoint->GetPosition();
G4ThreeVector hit_mom = preStepPoint->GetMomentumDirection();
double zint = hitPoint.z();
double pz = hit_mom.z();
// Check if theTrack moves backward
if (pz * zint < 0.) backward = true;
// Check that theTrack is above the energy threshold to use Shower Library
bool aboveThreshold = false;
if(theTrack->GetKineticEnergy() > energyThresholdSL) aboveThreshold = true;
// Check if theTrack is a muon (if so, DO NOT use Shower Library)
bool notaMuon = true;
G4int mumPDG = 13;
G4int mupPDG = -13;
G4int parCode = theTrack->GetDefinition()->GetPDGEncoding();
if (parCode == mupPDG || parCode == mumPDG ) notaMuon = false;
// angle condition
double theta_max = M_PI - 3.1305; // angle in radians corresponding to -5.2 eta
double R_mom=sqrt(hit_mom.x()*hit_mom.x() + hit_mom.y()*hit_mom.y());
double theta = atan2(R_mom,std::abs(pz));
bool angleok = false;
if ( theta < theta_max) angleok = true;
// OkToUse
double R = sqrt(hitPoint.x()*hitPoint.x() + hitPoint.y()*hitPoint.y());
bool dot = false;
if ( zint < -14450. && R < 45.) dot = true;
bool inRange = true;
if ( zint < -14700. || R > 193.) inRange = false;
bool OkToUse = false;
if ( inRange && !dot) OkToUse = true;
if (useShowerLibrary && aboveThreshold && notaMuon && (!backward) && OkToUse && angleok && currentLV == lvCAST ) {
// Use Castor shower library if energy is above threshold, is not a muon
// and is not moving backward
getFromLibrary(aStep);
#ifdef debugLog
LogDebug("ForwardSim") << " Current logical volume is " << nameVolume ;
#endif
return NCherPhot;
} else {
// Usual calculations
// G4ThreeVector hitPoint = preStepPoint->GetPosition();
// G4ThreeVector hit_mom = preStepPoint->GetMomentumDirection();
G4double stepl = aStep->GetStepLength()/cm;
G4double beta = preStepPoint->GetBeta();
G4double charge = preStepPoint->GetCharge();
// G4VProcess* curprocess = preStepPoint->GetProcessDefinedStep();
// G4String namePr = preStepPoint->GetProcessDefinedStep()->GetProcessName();
// std::string nameProcess;
// nameProcess.assign(namePr,0,4);
// G4LogicalVolume* lv = currentPV->GetLogicalVolume();
// G4Material* mat = lv->GetMaterial();
// G4double rad = mat->GetRadlen();
// postStepPoint information *********************************************
G4StepPoint* postStepPoint= aStep->GetPostStepPoint();
G4VPhysicalVolume* postPV = postStepPoint->GetPhysicalVolume();
G4String postname = postPV->GetName();
std::string postnameVolume;
postnameVolume.assign(postname,0,4);
// theTrack information *************************************************
// G4Track* theTrack = aStep->GetTrack();
//G4double entot = theTrack->GetTotalEnergy();
G4ThreeVector vert_mom = theTrack->GetVertexMomentumDirection();
G4ThreeVector localPoint = theTrack->GetTouchable()->GetHistory()->
GetTopTransform().TransformPoint(hitPoint);
// calculations... *************************************************
float phi = -100.;
if (vert_mom.x() != 0) phi = atan2(vert_mom.y(),vert_mom.x());
if (phi < 0.) phi += twopi;
G4String particleType = theTrack->GetDefinition()->GetParticleName();
#ifdef debugLog
float costheta =vert_mom.z()/sqrt(vert_mom.x()*vert_mom.x()+
vert_mom.y()*vert_mom.y()+
vert_mom.z()*vert_mom.z());
float theta = acos(std::min(std::max(costheta,float(-1.)),float(1.)));
float eta = -log(tan(theta/2));
G4int primaryID = theTrack->GetTrackID();
// *************************************************
// *************************************************
double edep = aStep->GetTotalEnergyDeposit();
#endif
// *************************************************
// *************************************************
// take into account light collection curve for plate
// double weight = curve_Castor(nameVolume, preStepPoint);
// double edep = aStep->GetTotalEnergyDeposit() * weight;
// *************************************************
// *************************************************
/* comments for sensitive volumes:
C001 ...-... CP06,CBU1 ...-...CALM --- > fibres and bundle
for first release of CASTOR
CASF --- > quartz plate for first and second releases of CASTOR
GF2Q, GFNQ, GR2Q, GRNQ
for tests with my own test geometry of HF (on ask of Gavrilov)
C3TF, C4TF - for third release of CASTOR
*/
double meanNCherPhot=0;
if (currentLV == lvC3EF || currentLV == lvC4EF || currentLV == lvC3HF ||
currentLV == lvC4HF) {
// if(nameVolume == "C3EF" || nameVolume == "C4EF" || nameVolume == "C3HF" || nameVolume == "C4HF") {
float bThreshold = 0.67;
float nMedium = 1.4925;
// float photEnSpectrDL = (1./400.nm-1./700.nm)*10000000.cm/nm; /* cm-1 */
// float photEnSpectrDL = 10714.285714;
float photEnSpectrDE = 1.24; /* see below */
/* E = 2pi*(1./137.)*(eV*cm/370.)/lambda = */
/* = 12.389184*(eV*cm)/lambda */
/* Emax = 12.389184*(eV*cm)/400nm*10-7cm/nm = 3.01 eV */
/* Emin = 12.389184*(eV*cm)/700nm*10-7cm/nm = 1.77 eV */
/* delE = Emax - Emin = 1.24 eV --> */
/* */
/* default for Castor nameVolume == "CASF" or (C3TF & C4TF) */
float thFullRefl = 23.; /* 23.dergee */
float thFullReflRad = thFullRefl*pi/180.;
/* default for Castor nameVolume == "CASF" or (C3TF & C4TF) */
float thFibDir = 45.; /* .dergee */
/* for test HF geometry volumes:
if(nameVolume == "GF2Q" || nameVolume == "GFNQ" ||
nameVolume == "GR2Q" || nameVolume == "GRNQ")
thFibDir = 0.0; // .dergee
*/
float thFibDirRad = thFibDir*pi/180.;
/* */
/* */
// at which theta the point is located:
// float th1 = hitPoint.theta();
// theta of charged particle in LabRF(hit momentum direction):
float costh =hit_mom.z()/sqrt(hit_mom.x()*hit_mom.x()+
hit_mom.y()*hit_mom.y()+
hit_mom.z()*hit_mom.z());
if (zint < 0) costh = -costh;
float th = acos(std::min(std::max(costh,float(-1.)),float(1.)));
// just in case (can do bot use):
if (th < 0.) th += twopi;
// theta of cone with Cherenkov photons w.r.t.direction of charged part.:
float costhcher =1./(nMedium*beta);
float thcher = acos(std::min(std::max(costhcher,float(-1.)),float(1.)));
// diff thetas of charged part. and quartz direction in LabRF:
float DelFibPart = fabs(th - thFibDirRad);
// define real distances:
float d = fabs(tan(th)-tan(thFibDirRad));
// float a = fabs(tan(thFibDirRad)-tan(thFibDirRad+thFullReflRad));
// float r = fabs(tan(th)-tan(th+thcher));
float a = tan(thFibDirRad)+tan(fabs(thFibDirRad-thFullReflRad));
float r = tan(th)+tan(fabs(th-thcher));
// define losses d_qz in cone of full reflection inside quartz direction
float d_qz;
#ifdef debugLog
float variant;
#endif
if(DelFibPart > (thFullReflRad + thcher) ) {
d_qz = 0.;
#ifdef debugLog
variant=0.;
#endif
}
// if(d > (r+a) ) {d_qz = 0.; variant=0.;}
else {
if((th + thcher) < (thFibDirRad+thFullReflRad) &&
(th - thcher) > (thFibDirRad-thFullReflRad)) {
d_qz = 1.;
#ifdef debugLog
variant=1.;
#endif
}
// if((DelFibPart + thcher) < thFullReflRad ) {d_qz = 1.; variant=1.;}
// if((d+r) < a ) {d_qz = 1.; variant=1.;}
else {
if((thFibDirRad + thFullReflRad) < (th + thcher) &&
(thFibDirRad - thFullReflRad) > (th - thcher) ) {
// if((thcher - DelFibPart ) > thFullReflRad )
// if((r-d ) > a )
d_qz = 0.;
#ifdef debugLog
variant=2.;
#endif
} else {
// if((thcher + DelFibPart ) > thFullReflRad &&
// thcher < (DelFibPart+thFullReflRad) )
// {
d_qz = 0.;
#ifdef debugLog
variant=3.;
#endif
// use crossed length of circles(cone projection)
// dC1/dC2 :
float arg_arcos = 0.;
float tan_arcos = 2.*a*d;
if(tan_arcos != 0.) arg_arcos =(r*r-a*a-d*d)/tan_arcos;
arg_arcos = fabs(arg_arcos);
float th_arcos = acos(std::min(std::max(arg_arcos,float(-1.)),float(1.)));
d_qz = th_arcos/pi/2.;
d_qz = fabs(d_qz);
// }
// else
// {
// d_qz = 0.; variant=4.;
//#ifdef debugLog
// std::cout <<" ===============>variant 4 information: <===== " <<std::endl;
// std::cout <<" !!!!!!!!!!!!!!!!!!!!!! variant = " << variant <<std::endl;
//#endif
//
// }
}
}
}
// !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
if(charge != 0. && beta > bThreshold ) {
meanNCherPhot = 370.*charge*charge*
( 1. - 1./(nMedium*nMedium*beta*beta) )*
photEnSpectrDE*stepl;
G4int poissNCherPhot = (G4int) G4Poisson(meanNCherPhot);
if(poissNCherPhot < 0) poissNCherPhot = 0;
float effPMTandTransport = 0.19;
double ReflPower = 0.1;
double proba = d_qz + (1-d_qz)*ReflPower;
NCherPhot = poissNCherPhot*effPMTandTransport*proba*0.307;
#ifdef debugLog
float thgrad = th*180./pi;
float thchergrad = thcher*180./pi;
float DelFibPartgrad = DelFibPart*180./pi;
LogDebug("ForwardSim") << " ==============================> start all "
<< "information:<========= \n" << " =====> for "
<< "test:<=== \n" << " variant = " << variant
<< "\n thgrad = " << thgrad << "\n thchergrad "
<< "= " << thchergrad << "\n DelFibPartgrad = "
<< DelFibPartgrad << "\n d_qz = " << d_qz
<< "\n =====> Start Step Information <=== \n"
<< " ===> calo preStepPoint info <=== \n"
<< " hitPoint = " << hitPoint << "\n"
<< " hitMom = " << hit_mom << "\n"
<< " stepControlFlag = " << stepControlFlag
// << "\n curprocess = " << curprocess << "\n"
// << " nameProcess = " << nameProcess
<< "\n charge = " << charge << "\n"
<< " beta = " << beta << "\n"
<< " bThreshold = " << bThreshold << "\n"
<< " thgrad =" << thgrad << "\n"
<< " effPMTandTransport=" << effPMTandTransport
// << "\n volume = " << name
<< "\n nameVolume = " << nameVolume << "\n"
<< " nMedium = " << nMedium << "\n"
// << " rad length = " << rad << "\n"
// << " material = " << mat << "\n"
<< " stepl = " << stepl << "\n"
<< " photEnSpectrDE = " << photEnSpectrDE <<"\n"
<< " edep = " << edep << "\n"
<< " ===> calo theTrack info <=== " << "\n"
<< " particleType = " << particleType << "\n"
<< " primaryID = " << primaryID << "\n"
<< " entot= " << theTrack->GetTotalEnergy() << "\n"
<< " vert_eta= " << eta << "\n"
<< " vert_phi= " << phi << "\n"
<< " vert_mom= " << vert_mom << "\n"
<< " ===> calo hit preStepPointinfo <=== "<<"\n"
<< " local point = " << localPoint << "\n"
<< " ==============================> final info"
<< ": <=== \n"
<< " meanNCherPhot = " << meanNCherPhot << "\n"
<< " poissNCherPhot = " << poissNCherPhot <<"\n"
<< " NCherPhot = " << NCherPhot;
#endif
// Included by WC
// std::cout << "\n volume = " << name
// << "\n nameVolume = " << nameVolume << "\n"
// << "\n postvolume = " << postname
// << "\n postnameVolume = " << postnameVolume << "\n"
// << "\n particleType = " << particleType
// << "\n primaryID = " << primaryID << "\n";
}
}
#ifdef debugLog
LogDebug("ForwardSim") << "CastorSD:: " << nameVolume
// << " Light Collection Efficiency " << weight
<< " Weighted Energy Deposit " << edep/MeV
<< " MeV\n";
#endif
// Temporary member for testing purpose only...
// unit_id = setDetUnitId(aStep);
// if(NCherPhot != 0) std::cout << "\n UnitID = " << unit_id << " ; NCherPhot = " << NCherPhot ;
return NCherPhot;
}
}
return 0;
}
| void CastorSD::getFromLibrary | ( | G4Step * | aStep | ) | [private] |
Definition at line 614 of file CastorSD.cc.
References CaloSD::checkHit(), CaloSD::createNewHit(), CaloSD::currentHit, CaloSD::currentID, CaloSD::edepositEM, CaloSD::edepositHAD, CaloSD::emPDG, CaloSD::epPDG, CaloSD::gammaPDG, CastorShowerEvent::getDetID(), CastorShowerEvent::getNhit(), CastorShowerEvent::getNphotons(), CastorShowerEvent::getPrimE(), CastorShowerLibrary::getShowerHits(), CastorShowerEvent::getTime(), GetVolume(), i, LogDebug, non_compensation_factor, convertSQLiteXML::ok, CaloSD::posGlobal, CaloSD::preStepPoint, CaloSD::previousID, CaloSD::resetForNewPrimary(), rotateUnitID(), pileupReCalc_HLTpaths::scale, CaloHitID::setID(), setTrackID(), showerLibrary, CaloSD::theTrack, cond::rpcobgas::time, funct::true, and CaloSD::updateHit().
Referenced by getEnergyDeposit().
{
//
// Method to get hits from the Shower Library
//
// CastorShowerEvent hits returned by getShowerHits are used to
// replace the full simulation of the shower from theTrack
//
// "updateHit" save the Hits to a CaloG4Hit container
//
preStepPoint = aStep->GetPreStepPoint();
theTrack = aStep->GetTrack();
bool ok;
// **** Call method to retrieve hits from the ShowerLibrary ****
CastorShowerEvent hits = showerLibrary->getShowerHits(aStep, ok);
double etrack = preStepPoint->GetKineticEnergy();
int primaryID = setTrackID(aStep);
// int primaryID = theTrack->GetTrackID();
// Reset entry point for new primary
posGlobal = preStepPoint->GetPosition();
resetForNewPrimary(posGlobal, etrack);
// Check whether track is EM or HAD
G4int particleCode = theTrack->GetDefinition()->GetPDGEncoding();
bool isEM , isHAD ;
if (particleCode==emPDG || particleCode==epPDG || particleCode==gammaPDG) {
isEM = true ; isHAD = false;
} else {
isEM = false; isHAD = true ;
}
#ifdef debugLog
// edm::LogInfo("ForwardSim") << "\n CastorSD::getFromLibrary: "
LogDebug("ForwardSim") << "\n CastorSD::getFromLibrary: "
<< hits.getNhit() << " hits for " << GetName() << " from "
<< theTrack->GetDefinition()->GetParticleName() << " of "
<< preStepPoint->GetKineticEnergy()/GeV << " GeV and trackID "
<< theTrack->GetTrackID() ;
#endif
// Scale to correct energy
double E_track = preStepPoint->GetTotalEnergy() ;
double E_SLhit = hits.getPrimE() * GeV ;
double scale = E_track/E_SLhit ;
//Non compensation
if (isHAD){
scale=scale*non_compensation_factor; // if hadronic extend the scale with the non-compensation factor
} else {
scale=scale; // if electromagnetic, don't do anything
}
/* double theTrackEnergy = theTrack->GetTotalEnergy() ;
if(fabs(theTrackEnergy-E_track)>10.) {
edm::LogInfo("ForwardSim") << "\n TrackID = " << theTrack->GetTrackID()
<< "\n theTrackEnergy = " << theTrackEnergy
<< "\n preStepPointEnergy = " << E_track ;
G4TrackVector tsec = *(aStep->GetSecondary());
for (unsigned int kk=0; kk<tsec.size(); kk++) {
edm::LogInfo("ForwardSim") << "CastorSD::getFromLibrary:"
<< "\n tsec[" << kk << "]->GetTrackID() = "
<< tsec[kk]->GetTrackID()
<< " with energy "
<< tsec[kk]->GetTotalEnergy() ;
}
}
*/
// Loop over hits retrieved from the library
for (unsigned int i=0; i<hits.getNhit(); i++) {
// Get nPhotoElectrons and set edepositEM / edepositHAD accordingly
double nPhotoElectrons = hits.getNphotons(i);
// Apply scaling
nPhotoElectrons *= scale ;
if(isEM) {
// edepositEM = nPhotoElectrons*GeV;
edepositEM = nPhotoElectrons;
edepositHAD = 0.;
} else if(isHAD) {
edepositEM = 0.;
edepositHAD = nPhotoElectrons;
// edepositHAD = nPhotoElectrons*GeV;
}
// Get hit position and time
double time = hits.getTime(i);
// math::XYZPoint position = hits.getHitPosition(i);
// Get hit detID
unsigned int unitID = hits.getDetID(i);
// Make the detID "rotation" from one sector to another taking into account the
// sectors of the impinging particle (theTrack) and of the particle that produced
// the 'hits' retrieved from shower library
unsigned int rotatedUnitID = rotateUnitID(unitID , theTrack , hits);
currentID.setID(rotatedUnitID, time, primaryID, 0);
// currentID.setID(unitID, time, primaryID, 0);
// check if it is in the same unit and timeslice as the previous one
if (currentID == previousID) {
updateHit(currentHit);
} else {
if (!checkHit()) currentHit = createNewHit();
}
} // End of loop over hits
//Now kill the current track
if (ok) {
theTrack->SetTrackStatus(fStopAndKill);
#ifdef debugLog
LogDebug("ForwardSim") << "CastorSD::getFromLibrary:"
<< "\n \"theTrack\" with TrackID() = "
<< theTrack->GetTrackID()
<< " and with energy "
<< theTrack->GetTotalEnergy()
<< " has been set to be killed" ;
#endif
G4TrackVector tv = *(aStep->GetSecondary());
for (unsigned int kk=0; kk<tv.size(); kk++) {
if (tv[kk]->GetVolume() == preStepPoint->GetPhysicalVolume()) {
tv[kk]->SetTrackStatus(fStopAndKill);
#ifdef debugLog
LogDebug("ForwardSim") << "CastorSD::getFromLibrary:"
<< "\n tv[" << kk << "]->GetTrackID() = "
<< tv[kk]->GetTrackID()
<< " with energy "
<< tv[kk]->GetTotalEnergy()
<< " has been set to be killed" ;
#endif
}
}
}
}
| void CastorSD::initRun | ( | ) | [protected, virtual] |
Reimplemented from CaloSD.
Definition at line 83 of file CastorSD.cc.
References CastorShowerLibrary::initParticleTable(), showerLibrary, and useShowerLibrary.
{
if (useShowerLibrary) {
// showerLibrary = new CastorShowerLibrary(name, cpv, p);
G4ParticleTable * theParticleTable = G4ParticleTable::GetParticleTable();
showerLibrary->initParticleTable(theParticleTable);
edm::LogInfo("ForwardSim") << "CastorSD::initRun: Using Castor Shower Library \n";
}
}
| uint32_t CastorSD::rotateUnitID | ( | uint32_t | unitID, |
| G4Track * | track, | ||
| CastorShowerEvent | shower | ||
| ) | [private] |
Definition at line 549 of file CastorSD.cc.
References printConversionInfo::aux, CastorShowerEvent::getPrimPhi(), LogDebug, M_PI, and reco::btau::trackPhi.
Referenced by getFromLibrary().
{
// ==============================================================
//
// o Exploit Castor phi symmetry to return newUnitID for
// shower hits based on track phi coordinate
//
// ==============================================================
// Get 'track' phi:
float trackPhi = track->GetPosition().phi();
if(trackPhi<0) trackPhi += 2*M_PI ;
// Get phi from primary that gave rise to SL 'shower':
float showerPhi = shower.getPrimPhi();
if(showerPhi<0) showerPhi += 2*M_PI ;
// Delta phi:
// Find the OctSector for which 'track' and 'shower' belong
int trackOctSector = (int) ( trackPhi / (M_PI/4) ) ;
int showerOctSector = (int) ( showerPhi / (M_PI/4) ) ;
uint32_t newUnitID;
uint32_t sec = ( ( unitID>>4 ) & 0xF ) ;
uint32_t complement = ( unitID & 0xFFFFFF0F ) ;
// edm::LogInfo("ForwardSim") << "\n CastorSD::rotateUnitID: "
// << "\n unitID = " << unitID
// << "\n sec = " << sec
// << "\n complement = " << complement ;
// Get 'track' z:
float trackZ = track->GetPosition().z();
int aux ;
int dSec = 2*(trackOctSector - showerOctSector) ;
// if(trackZ<0) // Good for revision 1.8 of CastorNumberingScheme
if(trackZ>0) // Good for revision 1.9 of CastorNumberingScheme
{
int sec1 = sec-dSec;
// sec -= dSec ;
if(sec1<0) sec1 += 16;
if(sec1>15) sec1 -= 16;
sec = (uint32_t)(sec1);
} else {
if( dSec<0 ) sec += 16 ;
sec += dSec ;
aux = (int) (sec/16) ;
sec -= aux*16 ;
}
sec = sec<<4 ;
newUnitID = complement | sec ;
#ifdef debugLog
if(newUnitID != unitID) {
LogDebug("ForwardSim") << "\n CastorSD::rotateUnitID: "
<< "\n unitID = " << unitID
<< "\n newUnitID = " << newUnitID ;
}
#endif
return newUnitID ;
}
| uint32_t CastorSD::setDetUnitId | ( | G4Step * | step | ) | [virtual] |
Implements CaloSD.
Definition at line 508 of file CastorSD.cc.
References CastorNumberingScheme::getUnitID(), and numberingScheme.
{
return (numberingScheme == 0 ? 0 : numberingScheme->getUnitID(aStep));
}
| void CastorSD::setNumberingScheme | ( | CastorNumberingScheme * | scheme | ) |
Definition at line 514 of file CastorSD.cc.
References numberingScheme.
Referenced by CastorSD().
{
if (scheme != 0) {
edm::LogInfo("ForwardSim") << "CastorSD: updates numbering scheme for "
<< GetName();
if (numberingScheme) delete numberingScheme;
numberingScheme = scheme;
}
}
| int CastorSD::setTrackID | ( | G4Step * | aStep | ) | [private] |
Definition at line 526 of file CastorSD.cc.
References TrackInformation::getIDonCaloSurface(), CaloSD::preStepPoint, CaloSD::previousID, CaloSD::resetForNewPrimary(), CaloSD::theTrack, and CaloHitID::trackID().
Referenced by getFromLibrary().
{
theTrack = aStep->GetTrack();
double etrack = preStepPoint->GetKineticEnergy();
TrackInformation * trkInfo = (TrackInformation *)(theTrack->GetUserInformation());
int primaryID = trkInfo->getIDonCaloSurface();
if (primaryID == 0) {
#ifdef debugLog
edm::LogWarning("ForwardSim") << "CastorSD: Problem with primaryID **** set by force "
<< "to TkID **** " << theTrack->GetTrackID();
#endif
primaryID = theTrack->GetTrackID();
}
if (primaryID != previousID.trackID())
resetForNewPrimary(preStepPoint->GetPosition(), etrack);
return primaryID;
}
double CastorSD::energyThresholdSL [private] |
Definition at line 52 of file CastorSD.h.
Referenced by CastorSD(), and getEnergyDeposit().
G4LogicalVolume* CastorSD::lvC3EF [private] |
Definition at line 48 of file CastorSD.h.
Referenced by CastorSD(), and getEnergyDeposit().
G4LogicalVolume * CastorSD::lvC3HF [private] |
Definition at line 48 of file CastorSD.h.
Referenced by CastorSD(), and getEnergyDeposit().
G4LogicalVolume * CastorSD::lvC4EF [private] |
Definition at line 48 of file CastorSD.h.
Referenced by CastorSD(), and getEnergyDeposit().
G4LogicalVolume * CastorSD::lvC4HF [private] |
Definition at line 48 of file CastorSD.h.
Referenced by CastorSD(), and getEnergyDeposit().
G4LogicalVolume* CastorSD::lvCAST [private] |
Definition at line 49 of file CastorSD.h.
Referenced by CastorSD(), and getEnergyDeposit().
double CastorSD::non_compensation_factor [private] |
Definition at line 53 of file CastorSD.h.
Referenced by CastorSD(), and getFromLibrary().
CastorNumberingScheme* CastorSD::numberingScheme [private] |
Definition at line 46 of file CastorSD.h.
Referenced by setDetUnitId(), and setNumberingScheme().
CastorShowerLibrary* CastorSD::showerLibrary [private] |
Definition at line 47 of file CastorSD.h.
Referenced by CastorSD(), getFromLibrary(), initRun(), and ~CastorSD().
bool CastorSD::useShowerLibrary [private] |
Definition at line 51 of file CastorSD.h.
Referenced by CastorSD(), getEnergyDeposit(), initRun(), and ~CastorSD().
1.7.3