File indexing completed on 2023-10-25 09:32:09
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "Alignment/LaserAlignmentSimulation/interface/LaserPrimaryGeneratorAction.h"
0010 #include "SimG4Core/Notification/interface/GenParticleInfo.h"
0011
0012 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0013
0014 #include "G4SystemOfUnits.hh"
0015
0016 LaserPrimaryGeneratorAction::LaserPrimaryGeneratorAction(edm::ParameterSet const &theConf)
0017 : thePhotonEnergy(0),
0018 thenParticleInGun(0),
0019 thenParticle(0),
0020 theLaserBeamsInTEC1(),
0021 theLaserBeamsInTEC2(),
0022 theLaserBeamsInTECTIBTOBTEC() {
0023
0024
0025
0026 thePhotonEnergy = theConf.getUntrackedParameter<double>("PhotonEnergy", 1.15) * eV;
0027
0028
0029 thenParticleInGun = theConf.getUntrackedParameter<int>("NumberOfPhotonsInParticleGun", 1);
0030
0031
0032
0033
0034
0035 thenParticle = theConf.getUntrackedParameter<int>("NumberOfPhotonsInEachBeam", 1);
0036
0037
0038
0039
0040
0041 theLaserBeamsInTEC1 = new LaserBeamsTEC1(thenParticleInGun, thenParticle, thePhotonEnergy);
0042
0043
0044 theLaserBeamsInTEC2 = new LaserBeamsTEC2(thenParticleInGun, thenParticle, thePhotonEnergy);
0045
0046
0047 theLaserBeamsInTECTIBTOBTEC = new LaserBeamsBarrel(thenParticleInGun, thenParticle, thePhotonEnergy);
0048
0049 }
0050
0051 LaserPrimaryGeneratorAction::~LaserPrimaryGeneratorAction() {
0052
0053
0054 if (theLaserBeamsInTEC1 != nullptr) {
0055 delete theLaserBeamsInTEC1;
0056 }
0057 if (theLaserBeamsInTEC2 != nullptr) {
0058 delete theLaserBeamsInTEC2;
0059 }
0060 if (theLaserBeamsInTECTIBTOBTEC != nullptr) {
0061 delete theLaserBeamsInTECTIBTOBTEC;
0062 }
0063
0064 }
0065
0066 void LaserPrimaryGeneratorAction::GeneratePrimaries(G4Event *myEvent) {
0067
0068
0069
0070
0071 LogDebug("LaserPrimaryGeneratorAction") << "<LaserPrimaryGeneratorAction::GeneratePrimaries(G4Event*)>: create a "
0072 "new Laser Event";
0073
0074
0075 theLaserBeamsInTEC1->GeneratePrimaries(myEvent);
0076
0077
0078 theLaserBeamsInTEC2->GeneratePrimaries(myEvent);
0079
0080
0081 theLaserBeamsInTECTIBTOBTEC->GeneratePrimaries(myEvent);
0082
0083
0084
0085 int theID = 0;
0086
0087 for (int i = 1; i < myEvent->GetNumberOfPrimaryVertex(); i++) {
0088 G4PrimaryVertex *theVertex = myEvent->GetPrimaryVertex(i);
0089
0090 for (int j = 0; j < theVertex->GetNumberOfParticle(); j++) {
0091 G4PrimaryParticle *thePrimary = theVertex->GetPrimary(j);
0092
0093 setGeneratorId(thePrimary, theID);
0094 theID++;
0095 }
0096 }
0097
0098 }
0099
0100 void LaserPrimaryGeneratorAction::setGeneratorId(G4PrimaryParticle *aParticle, int ID) const {
0101
0102
0103
0104
0105
0106
0107
0108 aParticle->SetUserInformation(new GenParticleInfo(ID));
0109
0110 }