From 464ce6f8957a3f07990a7ddd71784f3f54bd4201 Mon Sep 17 00:00:00 2001 From: Kamil Sedlak Date: Mon, 22 Jun 2009 16:39:50 +0000 Subject: [PATCH] 2009-06-22 Kamil Sedlak - introduced modifications done by Pavel Bakule at ISIS: a) the time when the initial muon is generated can be set according to a gaussian or flat distribution. (all muons were generated at time 0 before this change). b) a new volume "cylpart" has been added. The code can be compiled, but no tests of the newly implemented features have been done. In future, the information about the time of the muon generation should be added to the root output tree. --- include/musrPrimaryGeneratorAction.hh | 3 +++ include/musrPrimaryGeneratorMessenger.hh | 2 ++ src/musrDetectorConstruction.cc | 10 ++++++++++ src/musrPrimaryGeneratorAction.cc | 9 +++++++++ src/musrPrimaryGeneratorMessenger.cc | 16 ++++++++++++++++ 5 files changed, 40 insertions(+) diff --git a/include/musrPrimaryGeneratorAction.hh b/include/musrPrimaryGeneratorAction.hh index c5f00dd..d485d2c 100644 --- a/include/musrPrimaryGeneratorAction.hh +++ b/include/musrPrimaryGeneratorAction.hh @@ -52,6 +52,8 @@ class musrPrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction void Setvertex(G4ThreeVector v) {x0=v[0]; y0=v[1]; z0=v[2];} void SetvertexSigma(G4ThreeVector v) {xSigma=v[0]; ySigma=v[1]; zSigma=v[2];} void SetvertexBoundary(G4ThreeVector v) {rMaxAllowed=v[0]; zMinAllowed=v[1]; zMaxAllowed=v[2];} + void SetMuonTime(G4double val) {t0=val;} //P.B. 13 May 2009 + void SetMuonTimeSigma(G4double val) {tSigma=val;} //P.B. 13 May 2009 void SetKEnergy(G4double val); void SetMomentum(G4double val) {p0=val;} void SetMomentumSmearing(G4double val) {pSigma=val;} @@ -94,6 +96,7 @@ class musrPrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction static G4String thePrimaryParticleName ; G4double x0, y0, z0, xSigma, ySigma, zSigma, rMaxAllowed, zMinAllowed, zMaxAllowed; + G4double t0, tSigma; //P.B. 13 May 2009 G4double p0, pSigma, pMinAllowed, pMaxAllowed; G4double xangle0, yangle0, xangleSigma, yangleSigma, pitch; G4bool UnpolarisedMuonBeam, TransversalyUnpolarisedMuonBeam; diff --git a/include/musrPrimaryGeneratorMessenger.hh b/include/musrPrimaryGeneratorMessenger.hh index 858ffaf..a4f93be 100644 --- a/include/musrPrimaryGeneratorMessenger.hh +++ b/include/musrPrimaryGeneratorMessenger.hh @@ -51,6 +51,8 @@ class musrPrimaryGeneratorMessenger: public G4UImessenger G4UIcmdWith3VectorAndUnit* setvertexCmd; G4UIcmdWith3VectorAndUnit* setvertexSigmaCmd; G4UIcmdWith3VectorAndUnit* setvertexBoundaryCmd; + G4UIcmdWithADoubleAndUnit* setStarttimeCmd; //P.B. 13 May 2009 + G4UIcmdWithADoubleAndUnit* setStarttimeSigmaCmd; //P.B. 13 May 2009 G4UIcmdWithADoubleAndUnit* setKEnergyCmd; G4UIcmdWithADoubleAndUnit* setMomentumCmd; G4UIcmdWithADoubleAndUnit* setMomentumSmearingCmd; diff --git a/src/musrDetectorConstruction.cc b/src/musrDetectorConstruction.cc index 597e27a..617f9ec 100644 --- a/src/musrDetectorConstruction.cc +++ b/src/musrDetectorConstruction.cc @@ -216,6 +216,16 @@ G4VPhysicalVolume* musrDetectorConstruction::Construct() { solidName+=name; solid = new G4Para(solidName,x1*mm,x2*mm,x3*mm,x4*deg,x5*deg,x6*deg); } + else if (strcmp(tmpString2,"cylpart")==0){ // Volume introduced by Pavel Bakule on 12 May 2009 + sscanf(&line[0],"%*s %*s %*s %s %g %g %g %g %s %g %g %g %s %s", + name,&x1,&x2,&x3,&x4,material,&posx,&posy,&posz,mothersName,rotMatrix); + sscanf(&line[0],"%*s %*s %*s %*s %*g %*g %*g %*g %*s %*g %*g %*g %*s %*s %s %d %s",sensitiveDet,&volumeID,actualFieldName); + solidName+=name; + G4double roundingErr=0.01*mm; // to avoid some displaying problems of the subtracted volumes + G4Box* solidSubtractedBox = new G4Box("solidSubtractedBox",x2*mm,(x2-x4)*mm,x3*mm+roundingErr); + G4Tubs* solidCylinder = new G4Tubs("solidCylinder",0.*mm,x2*mm,x3*mm,0.*deg,180.*deg); + solid = new G4SubtractionSolid(solidName, solidCylinder, solidSubtractedBox); + } else if (strcmp(tmpString2,"uprofile")==0){ // Create a U-profile geometry. x1, x2, x3 define the outer dimensions of the U-profile (as a box), // x4 is the wall thickness of the U-profile. The centre of the U-profile diff --git a/src/musrPrimaryGeneratorAction.cc b/src/musrPrimaryGeneratorAction.cc index 2e5bb7a..542a6c3 100644 --- a/src/musrPrimaryGeneratorAction.cc +++ b/src/musrPrimaryGeneratorAction.cc @@ -58,6 +58,7 @@ musrPrimaryGeneratorAction::musrPrimaryGeneratorAction( musrDetectorConstruction* musrDC) :musrDetector(musrDC), x0(0), y0(0), z0(-10*cm), xSigma(0), ySigma(0), zSigma(0), rMaxAllowed(1e10*mm), zMinAllowed(-1e10*mm), zMaxAllowed(1e10*mm), + t0(0), tSigma(0), p0(0), pSigma(0), pMinAllowed(0), pMaxAllowed(1e10*mm), xangle0(0), yangle0(0), xangleSigma(0), yangleSigma(0), pitch(0), UnpolarisedMuonBeam(false), TransversalyUnpolarisedMuonBeam(false), xPolarisIni(1.), yPolarisIni(0.), zPolarisIni(0.), @@ -250,6 +251,13 @@ void musrPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) } // end of the part specific for the muons generated by random rather then from TURTLE + // Calculate particle (muon) starting time + G4double ParticleTime; //P.B. 13 May 2009 + if (tSigma>0) {ParticleTime = G4RandGauss::shoot(t0,tSigma);} // Gaussian distribution P.B. 13 May 2009 + else if (tSigma<0) {ParticleTime = t0 + tSigma*(G4UniformRand()*2.-1.);} // Uniform step distribution P.B. 13 May 2009 + else {ParticleTime = t0;} // Point-like P.B. 13 May 2009 + + // Calculate the final momentum G4double px, py, pz; px = p*sin(xangle); @@ -284,6 +292,7 @@ void musrPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) G4double particle_mass = particleGun->GetParticleDefinition()->GetPDGMass(); G4double particleEnergy = std::sqrt(p*p+particle_mass*particle_mass)-particle_mass; particleGun->SetParticleEnergy(particleEnergy); + particleGun->SetParticleTime(ParticleTime); //P.B. 13 May 2009 particleGun->SetParticleMomentumDirection(G4ThreeVector(px,py,pz)); particleGun->SetParticlePolarization(G4ThreeVector(xpolaris,ypolaris,zpolaris)); particleGun->GeneratePrimaryVertex(anEvent); diff --git a/src/musrPrimaryGeneratorMessenger.cc b/src/musrPrimaryGeneratorMessenger.cc index e6bb5f0..8bcc112 100644 --- a/src/musrPrimaryGeneratorMessenger.cc +++ b/src/musrPrimaryGeneratorMessenger.cc @@ -54,6 +54,16 @@ musrPrimaryGeneratorMessenger::musrPrimaryGeneratorMessenger(musrPrimaryGenerato setvertexBoundaryCmd->SetParameterName("mes_rMaxAllowed","mes_zMinAllowed","mes_zMaxAllowed",true,true); setvertexBoundaryCmd->SetDefaultUnit("mm"); + setStarttimeCmd = new G4UIcmdWithADoubleAndUnit("/gun/starttime",this); //P.B. 13 May 2009 + setStarttimeCmd->SetGuidance(" Set start time t of the generated muons (with unit)"); //P.B. 13 May 2009 + setStarttimeCmd->SetParameterName("mes_t0",true); //P.B. 13 May 2009 + setStarttimeCmd->SetDefaultUnit("ns"); //P.B. 13 May 2009 + + setStarttimeSigmaCmd = new G4UIcmdWithADoubleAndUnit("/gun/starttimesigma",this); //P.B. 13 May 2009 + setStarttimeSigmaCmd->SetGuidance(" Set start time sigma tSigma of the generated muons (with unit)"); //P.B. 13 May 2009 + setStarttimeSigmaCmd->SetParameterName("mes_tSigma",true); //P.B. 13 May 2009 + setStarttimeSigmaCmd->SetDefaultUnit("ns"); //P.B. 13 May 2009 + setKEnergyCmd = new G4UIcmdWithADoubleAndUnit("/gun/kenergy",this); setKEnergyCmd->SetGuidance(" Set kinetic energy of the generated muons (with unit)"); setKEnergyCmd->SetParameterName("mes_E0",true); @@ -144,6 +154,8 @@ musrPrimaryGeneratorMessenger::~musrPrimaryGeneratorMessenger() delete setvertexCmd; delete setvertexSigmaCmd; delete setvertexBoundaryCmd; + delete setStarttimeCmd; //P.B. 13 May 2009 + delete setStarttimeSigmaCmd; //P.B. 13 May 2009 delete setKEnergyCmd; delete setMomentumCmd; delete setMomentumSmearingCmd; @@ -172,6 +184,10 @@ void musrPrimaryGeneratorMessenger::SetNewValue(G4UIcommand * command,G4String n { musrAction->SetvertexSigma(setvertexSigmaCmd->GetNew3VectorValue(newValue));} if( command == setvertexBoundaryCmd) { musrAction->SetvertexBoundary(setvertexBoundaryCmd->GetNew3VectorValue(newValue));} + if( command == setStarttimeCmd) //P.B. 13 May 2009 + { musrAction->SetMuonTime(setStarttimeCmd->GetNewDoubleValue(newValue));} //P.B. 13 May 2009 + if( command == setStarttimeSigmaCmd) //P.B. 13 May 2009 + { musrAction->SetMuonTimeSigma(setStarttimeSigmaCmd->GetNewDoubleValue(newValue));} //P.B. 13 May 2009 if( command == setKEnergyCmd) { musrAction->SetKEnergy(setKEnergyCmd->GetNewDoubleValue(newValue));} if( command == setMomentumCmd)