#include <MaterialEffects.h>
Public Member Functions | |
double | energyLoss () const |
Return the energy loss by ionization in the current layer. | |
EnergyLossSimulator * | energyLossSimulator () const |
Return the Energy Loss engine. | |
void | interact (FSimEvent &simEvent, const TrackerLayer &layer, ParticlePropagator &PP, unsigned i) |
MaterialEffects (const edm::ParameterSet &matEff, const RandomEngine *engine) | |
Constructor. | |
MultipleScatteringSimulator * | multipleScatteringSimulator () const |
Return the Multiple Scattering engine. | |
MuonBremsstrahlungSimulator * | muonBremsstrahlungSimulator () const |
Return the Muon Bremsstrahlung engine. | |
void | save () |
Save nuclear interaction information. | |
double | thickness () const |
Return the thickness of the current layer. | |
~MaterialEffects () | |
Default destructor. | |
Private Member Functions | |
GlobalVector | normalVector (const TrackerLayer &layer, ParticlePropagator &myTrack) const |
The vector normal to the surface traversed. | |
double | radLengths (const TrackerLayer &layer, ParticlePropagator &myTrack) |
The number of radiation lengths traversed. | |
Private Attributes | |
BremsstrahlungSimulator * | Bremsstrahlung |
EnergyLossSimulator * | EnergyLoss |
MultipleScatteringSimulator * | MultipleScattering |
MuonBremsstrahlungSimulator * | MuonBremsstrahlung |
NuclearInteractionSimulator * | NuclearInteraction |
PairProductionSimulator * | PairProduction |
double | pTmin |
const RandomEngine * | random |
double | theEnergyLoss |
GlobalVector | theNormalVector |
double | theTECFudgeFactor |
double | theThickness |
bool | use_hardcoded |
Definition at line 50 of file MaterialEffects.h.
MaterialEffects::MaterialEffects | ( | const edm::ParameterSet & | matEff, |
const RandomEngine * | engine | ||
) |
Constructor.
Definition at line 25 of file MaterialEffects.cc.
References funct::A, Bremsstrahlung, EnergyLoss, edm::ParameterSet::getParameter(), edm::ParameterSet::getUntrackedParameter(), i, analyzePatCOC_cfg::inputFile, j, MultipleScattering, MuonBremsstrahlung, NuclearInteraction, PairProduction, pTmin, random, plotFactory::ratios, AlCaHLTBitMon_QueryRunRegistry::string, theTECFudgeFactor, use_hardcoded, and Gflash::Z.
: PairProduction(0), Bremsstrahlung(0),MuonBremsstrahlung(0), MultipleScattering(0), EnergyLoss(0), NuclearInteraction(0), pTmin(999.), random(engine), use_hardcoded(1) { // Set the minimal photon energy for a Brem from e+/- use_hardcoded = matEff.getParameter<bool >("use_hardcoded_geometry"); bool doPairProduction = matEff.getParameter<bool>("PairProduction"); bool doBremsstrahlung = matEff.getParameter<bool>("Bremsstrahlung"); bool doEnergyLoss = matEff.getParameter<bool>("EnergyLoss"); bool doMultipleScattering = matEff.getParameter<bool>("MultipleScattering"); bool doNuclearInteraction = matEff.getParameter<bool>("NuclearInteraction"); bool doMuonBremsstrahlung = matEff.getParameter<bool>("MuonBremsstrahlung"); double A = matEff.getParameter<double>("A"); double Z = matEff.getParameter<double>("Z"); double density = matEff.getParameter<double>("Density"); double radLen = matEff.getParameter<double>("RadiationLength"); // Set the minimal pT before giving up the dE/dx treatment if ( doPairProduction ) { double photonEnergy = matEff.getParameter<double>("photonEnergy"); PairProduction = new PairProductionSimulator(photonEnergy, random); } if ( doBremsstrahlung ) { double bremEnergy = matEff.getParameter<double>("bremEnergy"); double bremEnergyFraction = matEff.getParameter<double>("bremEnergyFraction"); Bremsstrahlung = new BremsstrahlungSimulator(bremEnergy, bremEnergyFraction, random); } //muon Brem+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ if ( doMuonBremsstrahlung ) { double bremEnergy = matEff.getParameter<double>("bremEnergy"); double bremEnergyFraction = matEff.getParameter<double>("bremEnergyFraction"); MuonBremsstrahlung = new MuonBremsstrahlungSimulator(random,A,Z,density,radLen,bremEnergy, bremEnergyFraction); } //++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ if ( doEnergyLoss ) { pTmin = matEff.getParameter<double>("pTmin"); EnergyLoss = new EnergyLossSimulator(random,A,Z,density,radLen); } if ( doMultipleScattering ) { MultipleScattering = new MultipleScatteringSimulator(random,A,Z,density,radLen); } if ( doNuclearInteraction ) { // The energies simulated std::vector<double> pionEnergies = matEff.getUntrackedParameter<std::vector<double> >("pionEnergies"); // The particle types simulated std::vector<int> pionTypes = matEff.getUntrackedParameter<std::vector<int> >("pionTypes"); // The corresponding particle names std::vector<std::string> pionNames = matEff.getUntrackedParameter<std::vector<std::string> >("pionNames"); // The corresponding particle masses std::vector<double> pionMasses = matEff.getUntrackedParameter<std::vector<double> >("pionMasses"); // The smallest momentum for inelastic interactions std::vector<double> pionPMin = matEff.getUntrackedParameter<std::vector<double> >("pionMinP"); // The interaction length / radiation length ratio for each particle type std::vector<double> lengthRatio = matEff.getParameter<std::vector<double> >("lengthRatio"); // std::map<int,double> lengthRatio; // for ( unsigned i=0; i<theLengthRatio.size(); ++i ) // lengthRatio[ pionTypes[i] ] = theLengthRatio[i]; // A global fudge factor for TEC layers (which apparently do not react to // hadrons the same way as all other layers... theTECFudgeFactor = matEff.getParameter<double>("fudgeFactor"); // The evolution of the interaction lengths with energy std::vector<double> theRatios = matEff.getUntrackedParameter<std::vector<double> >("ratios"); //std::map<int,std::vector<double> > ratios; //for ( unsigned i=0; i<pionTypes.size(); ++i ) { // for ( unsigned j=0; j<pionEnergies.size(); ++j ) { // ratios[ pionTypes[i] ].push_back(theRatios[ i*pionEnergies.size() + j ]); // } //} std::vector< std::vector<double> > ratios; ratios.resize(pionTypes.size()); for ( unsigned i=0; i<pionTypes.size(); ++i ) { for ( unsigned j=0; j<pionEnergies.size(); ++j ) { ratios[i].push_back(theRatios[ i*pionEnergies.size() + j ]); } } // The smallest momentum for elastic interactions double pionEnergy = matEff.getParameter<double>("pionEnergy"); // The algorithm to compute the distance between primary and secondaries // when a nuclear interaction occurs unsigned distAlgo = matEff.getParameter<unsigned>("distAlgo"); double distCut = matEff.getParameter<double>("distCut"); // The file to read the starting interaction in each files // (random reproducibility in case of a crash) std::string inputFile = matEff.getUntrackedParameter<std::string>("inputFile"); // Build the ID map (i.e., what is to be considered as a proton, etc...) std::map<int,int> idMap; // Protons std::vector<int> idProtons = matEff.getUntrackedParameter<std::vector<int> >("protons"); for ( unsigned i=0; i<idProtons.size(); ++i ) idMap[idProtons[i]] = 2212; // Anti-Protons std::vector<int> idAntiProtons = matEff.getUntrackedParameter<std::vector<int> >("antiprotons"); for ( unsigned i=0; i<idAntiProtons.size(); ++i ) idMap[idAntiProtons[i]] = -2212; // Neutrons std::vector<int> idNeutrons = matEff.getUntrackedParameter<std::vector<int> >("neutrons"); for ( unsigned i=0; i<idNeutrons.size(); ++i ) idMap[idNeutrons[i]] = 2112; // Anti-Neutrons std::vector<int> idAntiNeutrons = matEff.getUntrackedParameter<std::vector<int> >("antineutrons"); for ( unsigned i=0; i<idAntiNeutrons.size(); ++i ) idMap[idAntiNeutrons[i]] = -2112; // K0L's std::vector<int> idK0Ls = matEff.getUntrackedParameter<std::vector<int> >("K0Ls"); for ( unsigned i=0; i<idK0Ls.size(); ++i ) idMap[idK0Ls[i]] = 130; // K+'s std::vector<int> idKplusses = matEff.getUntrackedParameter<std::vector<int> >("Kplusses"); for ( unsigned i=0; i<idKplusses.size(); ++i ) idMap[idKplusses[i]] = 321; // K-'s std::vector<int> idKminusses = matEff.getUntrackedParameter<std::vector<int> >("Kminusses"); for ( unsigned i=0; i<idKminusses.size(); ++i ) idMap[idKminusses[i]] = -321; // pi+'s std::vector<int> idPiplusses = matEff.getUntrackedParameter<std::vector<int> >("Piplusses"); for ( unsigned i=0; i<idPiplusses.size(); ++i ) idMap[idPiplusses[i]] = 211; // pi-'s std::vector<int> idPiminusses = matEff.getUntrackedParameter<std::vector<int> >("Piminusses"); for ( unsigned i=0; i<idPiminusses.size(); ++i ) idMap[idPiminusses[i]] = -211; // Construction NuclearInteraction = new NuclearInteractionSimulator(pionEnergies, pionTypes, pionNames, pionMasses, pionPMin, pionEnergy, lengthRatio, ratios, idMap, inputFile, distAlgo, distCut, random); } }
MaterialEffects::~MaterialEffects | ( | ) |
Default destructor.
Definition at line 218 of file MaterialEffects.cc.
References Bremsstrahlung, EnergyLoss, MultipleScattering, MuonBremsstrahlung, NuclearInteraction, and PairProduction.
{ if ( PairProduction ) delete PairProduction; if ( Bremsstrahlung ) delete Bremsstrahlung; if ( EnergyLoss ) delete EnergyLoss; if ( MultipleScattering ) delete MultipleScattering; if ( NuclearInteraction ) delete NuclearInteraction; //Muon Brem if ( MuonBremsstrahlung ) delete MuonBremsstrahlung; }
double MaterialEffects::energyLoss | ( | ) | const [inline] |
Return the energy loss by ionization in the current layer.
Definition at line 76 of file MaterialEffects.h.
References theEnergyLoss.
Referenced by TrajectoryManager::createPSimHits().
{ return theEnergyLoss; }
EnergyLossSimulator* MaterialEffects::energyLossSimulator | ( | ) | const [inline] |
Return the Energy Loss engine.
Definition at line 84 of file MaterialEffects.h.
References EnergyLoss.
Referenced by MuonSimHitProducer::applyMaterialEffects(), and CalorimetryManager::MuonMipSimulation().
{ return EnergyLoss; }
void MaterialEffects::interact | ( | FSimEvent & | simEvent, |
const TrackerLayer & | layer, | ||
ParticlePropagator & | PP, | ||
unsigned | i | ||
) |
Steer the various interaction processes in the Tracker Material and update the FSimEvent
Energy loss ---------------
Multiple scattering -----------------------
Definition at line 229 of file MaterialEffects.cc.
References abs, FBaseSimEvent::addSimTrack(), FBaseSimEvent::addSimVertex(), MaterialEffectsSimulator::beginDaughters(), FSimVertexType::BREM_VERTEX, Bremsstrahlung, RawParticle::charge(), MaterialEffectsSimulator::closestDaughterId(), MaterialEffectsSimulator::endDaughters(), EnergyLoss, TrackerLayer::layerNumber(), MultipleScattering, MuonBremsstrahlung, MaterialEffectsSimulator::nDaughters(), normalVector(), FSimVertexType::NUCL_VERTEX, NuclearInteraction, FSimVertexType::PAIR_VERTEX, PairProduction, RawParticle::pid(), FSimVertex::position(), pTmin, radLengths(), FSimTrack::setClosestDaughterId(), MaterialEffectsSimulator::setNormalVector(), theEnergyLoss, theNormalVector, theTECFudgeFactor, FBaseSimEvent::track(), MaterialEffectsSimulator::updateState(), use_hardcoded, RawParticle::vertex(), and FSimTrack::vertex().
Referenced by TrajectoryManager::reconstruct().
{ MaterialEffectsSimulator::RHEP_const_iter DaughterIter; double radlen; theEnergyLoss = 0; theNormalVector = normalVector(layer,myTrack); radlen = radLengths(layer,myTrack); //------------------- // Photon Conversion //------------------- if ( PairProduction && myTrack.pid()==22 ) { // PairProduction->updateState(myTrack,radlen); if ( PairProduction->nDaughters() ) { //add a vertex to the mother particle int ivertex = mySimEvent.addSimVertex(myTrack.vertex(),itrack, FSimVertexType::PAIR_VERTEX); // Check if it is a valid vertex first: if (ivertex>=0) { // This was a photon that converted for ( DaughterIter = PairProduction->beginDaughters(); DaughterIter != PairProduction->endDaughters(); ++DaughterIter) { mySimEvent.addSimTrack(&(*DaughterIter), ivertex); } // The photon converted. Return. return; } else { edm::LogWarning("MaterialEffects") << " WARNING: A non valid vertex was found in photon conv. -> " << ivertex << std::endl; } } } if ( myTrack.pid() == 22 ) return; //------------------------ // Nuclear interactions //------------------------ if ( NuclearInteraction && abs(myTrack.pid()) > 100 && abs(myTrack.pid()) < 1000000) { // Simulate a nuclear interaction double factor = 1.0; if(use_hardcoded){ if (layer.layerNumber() >= 19 && layer.layerNumber() <= 27 ) factor = theTECFudgeFactor; } NuclearInteraction->updateState(myTrack,radlen*factor); if ( NuclearInteraction->nDaughters() ) { //add a end vertex to the mother particle int ivertex = mySimEvent.addSimVertex(myTrack.vertex(),itrack, FSimVertexType::NUCL_VERTEX); // Check if it is a valid vertex first: if (ivertex>=0) { // This was a hadron that interacted inelastically int idaugh = 0; for ( DaughterIter = NuclearInteraction->beginDaughters(); DaughterIter != NuclearInteraction->endDaughters(); ++DaughterIter) { // The daughter in the event int daughId = mySimEvent.addSimTrack(&(*DaughterIter), ivertex); // Store the closest daughter in the mother info (for later tracking purposes) if ( NuclearInteraction->closestDaughterId() == idaugh++ ) { if ( mySimEvent.track(itrack).vertex().position().Pt() < 4.0 ) mySimEvent.track(itrack).setClosestDaughterId(daughId); } } // The hadron is destroyed. Return. return; } else { edm::LogWarning("MaterialEffects") << " WARNING: A non valid vertex was found in nucl. int. -> " << ivertex << std::endl; } } } if ( myTrack.charge() == 0 ) return; if ( !MuonBremsstrahlung && !Bremsstrahlung && !EnergyLoss && !MultipleScattering ) return; //---------------- // Bremsstrahlung //---------------- if ( Bremsstrahlung && abs(myTrack.pid())==11 ) { Bremsstrahlung->updateState(myTrack,radlen); if ( Bremsstrahlung->nDaughters() ) { // Add a vertex, but do not attach it to the electron, because it // continues its way... int ivertex = mySimEvent.addSimVertex(myTrack.vertex(),itrack, FSimVertexType::BREM_VERTEX); // Check if it is a valid vertex first: if (ivertex>=0) { for ( DaughterIter = Bremsstrahlung->beginDaughters(); DaughterIter != Bremsstrahlung->endDaughters(); ++DaughterIter) { mySimEvent.addSimTrack(&(*DaughterIter), ivertex); } } else { edm::LogWarning("MaterialEffects") << " WARNING: A non valid vertex was found in brem -> " << ivertex << std::endl; } } } //--------------------------- // Muon_Bremsstrahlung //-------------------------- if ( MuonBremsstrahlung && abs(myTrack.pid())==13 ) { MuonBremsstrahlung->updateState(myTrack,radlen); if ( MuonBremsstrahlung->nDaughters() ) { // Add a vertex, but do not attach it to the muon, because it // continues its way... int ivertex = mySimEvent.addSimVertex(myTrack.vertex(),itrack, FSimVertexType::BREM_VERTEX); // Check if it is a valid vertex first: if (ivertex>=0) { for ( DaughterIter = MuonBremsstrahlung->beginDaughters(); DaughterIter != MuonBremsstrahlung->endDaughters(); ++DaughterIter) { mySimEvent.addSimTrack(&(*DaughterIter), ivertex); } } else { edm::LogWarning("MaterialEffects") << " WARNING: A non valid vertex was found in muon brem -> " << ivertex << std::endl; } } } if ( EnergyLoss ) { theEnergyLoss = myTrack.E(); EnergyLoss->updateState(myTrack,radlen); theEnergyLoss -= myTrack.E(); } if ( MultipleScattering && myTrack.Pt() > pTmin ) { // MultipleScattering->setNormalVector(normalVector(layer,myTrack)); MultipleScattering->setNormalVector(theNormalVector); MultipleScattering->updateState(myTrack,radlen); } }
MultipleScatteringSimulator* MaterialEffects::multipleScatteringSimulator | ( | ) | const [inline] |
Return the Multiple Scattering engine.
Definition at line 79 of file MaterialEffects.h.
References MultipleScattering.
Referenced by MuonSimHitProducer::applyMaterialEffects().
{ return MultipleScattering; }
MuonBremsstrahlungSimulator* MaterialEffects::muonBremsstrahlungSimulator | ( | ) | const [inline] |
Return the Muon Bremsstrahlung engine.
Definition at line 89 of file MaterialEffects.h.
References MuonBremsstrahlung.
Referenced by MuonSimHitProducer::applyMaterialEffects().
{ return MuonBremsstrahlung; }
GlobalVector MaterialEffects::normalVector | ( | const TrackerLayer & | layer, |
ParticlePropagator & | myTrack | ||
) | const [private] |
The vector normal to the surface traversed.
Definition at line 459 of file MaterialEffects.cc.
References TrackerLayer::disk(), TrackerLayer::forward(), RawParticle::R(), RawParticle::X(), and RawParticle::Y().
Referenced by interact().
double MaterialEffects::radLengths | ( | const TrackerLayer & | layer, |
ParticlePropagator & | myTrack | ||
) | [private] |
The number of radiation lengths traversed.
Definition at line 418 of file MaterialEffects.cc.
References TrackerLayer::forward(), TrackerLayer::fudgeFactor(), TrackerLayer::fudgeMax(), TrackerLayer::fudgeMin(), TrackerLayer::fudgeNumber(), P, RawParticle::R(), TrackerLayer::surface(), theNormalVector, theThickness, and RawParticle::Z().
Referenced by interact().
{ // Thickness of layer theThickness = layer.surface().mediumProperties().radLen(); GlobalVector P(myTrack.Px(),myTrack.Py(),myTrack.Pz()); // Effective length of track inside layer (considering crossing angle) // double radlen = theThickness / fabs(P.dot(theNormalVector)/(P.mag()*theNormalVector.mag())); double radlen = theThickness / fabs(P.dot(theNormalVector)) * P.mag(); // This is a series of fudge factors (from the geometry description), // to describe the layer inhomogeneities (services, cables, supports...) double rad = myTrack.R(); double zed = fabs(myTrack.Z()); double factor = 1; // Are there fudge factors for this layer if ( layer.fudgeNumber() ) // If yes, loop on them for ( unsigned int iLayer=0; iLayer < layer.fudgeNumber(); ++iLayer ) { // Apply to R if forward layer, to Z if barrel layer if ( ( layer.forward() && layer.fudgeMin(iLayer) < rad && rad < layer.fudgeMax(iLayer) ) || ( !layer.forward() && layer.fudgeMin(iLayer) < zed && zed < layer.fudgeMax(iLayer) ) ) { factor = layer.fudgeFactor(iLayer); break; } } theThickness *= factor; return radlen * factor; }
void MaterialEffects::save | ( | ) |
Save nuclear interaction information.
Definition at line 467 of file MaterialEffects.cc.
References NuclearInteraction, and NuclearInteractionSimulator::save().
Referenced by TrajectoryManager::reconstruct().
{ // Save current nuclear interactions in the event libraries. if ( NuclearInteraction ) NuclearInteraction->save(); }
double MaterialEffects::thickness | ( | ) | const [inline] |
Return the thickness of the current layer.
Definition at line 73 of file MaterialEffects.h.
References theThickness.
Referenced by TrajectoryManager::createPSimHits().
{ return theThickness; }
Definition at line 106 of file MaterialEffects.h.
Referenced by interact(), MaterialEffects(), and ~MaterialEffects().
EnergyLossSimulator* MaterialEffects::EnergyLoss [private] |
Definition at line 110 of file MaterialEffects.h.
Referenced by energyLossSimulator(), interact(), MaterialEffects(), and ~MaterialEffects().
Definition at line 109 of file MaterialEffects.h.
Referenced by interact(), MaterialEffects(), multipleScatteringSimulator(), and ~MaterialEffects().
Definition at line 108 of file MaterialEffects.h.
Referenced by interact(), MaterialEffects(), muonBremsstrahlungSimulator(), and ~MaterialEffects().
Definition at line 111 of file MaterialEffects.h.
Referenced by interact(), MaterialEffects(), save(), and ~MaterialEffects().
Definition at line 105 of file MaterialEffects.h.
Referenced by interact(), MaterialEffects(), and ~MaterialEffects().
double MaterialEffects::pTmin [private] |
Definition at line 114 of file MaterialEffects.h.
Referenced by interact(), and MaterialEffects().
const RandomEngine* MaterialEffects::random [private] |
Definition at line 124 of file MaterialEffects.h.
Referenced by MaterialEffects().
double MaterialEffects::theEnergyLoss [private] |
Definition at line 117 of file MaterialEffects.h.
Referenced by energyLoss(), and interact().
GlobalVector MaterialEffects::theNormalVector [private] |
Definition at line 115 of file MaterialEffects.h.
Referenced by interact(), and radLengths().
double MaterialEffects::theTECFudgeFactor [private] |
Definition at line 118 of file MaterialEffects.h.
Referenced by interact(), and MaterialEffects().
double MaterialEffects::theThickness [private] |
Definition at line 116 of file MaterialEffects.h.
Referenced by radLengths(), and thickness().
bool MaterialEffects::use_hardcoded [private] |
Definition at line 126 of file MaterialEffects.h.
Referenced by interact(), and MaterialEffects().