LMU/musrsim
5
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

musrsim version which compiles/links with Geant4.10

Browse Source
This commit is contained in:
nemu
2014-08-15 16:34:18 +02:00
parent 2a200dedba
commit f60f4874d6
4 changed files with 18 additions and 17 deletions
Download Patch File Download Diff File Expand all files Collapse all files

1
include/MuDecayChannel.hh
View File

@ -38,6 +38,7 @@
#ifndef MuDecayChannel_h #ifndef MuDecayChannel_h
#define MuDecayChannel_h 1 #define MuDecayChannel_h 1
#include <CLHEP/Units/PhysicalConstants.h>
#include "G4ios.hh" #include "G4ios.hh"
#include "globals.hh" #include "globals.hh"
#include "G4VDecayChannel.hh" #include "G4VDecayChannel.hh"

1
include/musrRootOutput.hh
View File

@ -24,6 +24,7 @@
#ifndef musrRootOutput_h #ifndef musrRootOutput_h
#define musrRootOutput_h 1 #define musrRootOutput_h 1
//#include "G4UserRunAction.hh" //#include "G4UserRunAction.hh"
#include <CLHEP/Units/PhysicalConstants.h>
#include "globals.hh" #include "globals.hh"
#include "G4ThreeVector.hh" #include "G4ThreeVector.hh"
// ROOT // ROOT

21
src/MuDecayChannelWithSpin.cc
View File

@ -45,7 +45,6 @@
// ------------------------------------------------------------ // ------------------------------------------------------------
// //
#include "MuDecayChannelWithSpin.hh" #include "MuDecayChannelWithSpin.hh"
#include "Randomize.hh" #include "Randomize.hh"
#include "G4DecayProducts.hh" #include "G4DecayProducts.hh"
#include "G4LorentzVector.hh" #include "G4LorentzVector.hh"
@ -70,11 +69,11 @@ G4DecayProducts *MuDecayChannelWithSpin::DecayIt(G4double)
if (GetVerboseLevel()>1) G4cout << "MuDecayChannelWithSpin::DecayIt "; if (GetVerboseLevel()>1) G4cout << "MuDecayChannelWithSpin::DecayIt ";
#endif #endif
if (parent == 0) FillParent(); if (G4MT_parent == 0) FillParent();
if (daughters == 0) FillDaughters(); if (G4MT_daughters == 0) FillDaughters();
// parent mass // parent mass
G4double parentmass = parent->GetPDGMass(); G4double parentmass = G4MT_parent->GetPDGMass();
EMMU = parentmass; EMMU = parentmass;
@ -82,7 +81,7 @@ G4DecayProducts *MuDecayChannelWithSpin::DecayIt(G4double)
G4double daughtermass[3]; G4double daughtermass[3];
G4double sumofdaughtermass = 0.0; G4double sumofdaughtermass = 0.0;
for (G4int index=0; index<3; index++){ for (G4int index=0; index<3; index++){
daughtermass[index] = daughters[index]->GetPDGMass(); daughtermass[index] = G4MT_daughters[index]->GetPDGMass();
sumofdaughtermass += daughtermass[index]; sumofdaughtermass += daughtermass[index];
} }
@ -90,7 +89,7 @@ G4DecayProducts *MuDecayChannelWithSpin::DecayIt(G4double)
//create parent G4DynamicParticle at rest //create parent G4DynamicParticle at rest
G4ThreeVector dummy; G4ThreeVector dummy;
G4DynamicParticle * parentparticle = new G4DynamicParticle( parent, dummy, 0.0); G4DynamicParticle * parentparticle = new G4DynamicParticle( G4MT_parent, dummy, 0.0);
//create G4Decayproducts //create G4Decayproducts
G4DecayProducts *products = new G4DecayProducts(*parentparticle); G4DecayProducts *products = new G4DecayProducts(*parentparticle);
delete parentparticle; delete parentparticle;
@ -176,7 +175,7 @@ G4DecayProducts *MuDecayChannelWithSpin::DecayIt(G4double)
rndm = G4UniformRand(); rndm = G4UniformRand();
G4double phi = twopi * rndm; G4double phi = CLHEP::twopi * rndm;
if(energy < EMASS) energy = EMASS; if(energy < EMASS) energy = EMASS;
@ -200,7 +199,7 @@ G4DecayProducts *MuDecayChannelWithSpin::DecayIt(G4double)
direction0.rotateUz(parent_polarization); direction0.rotateUz(parent_polarization);
G4DynamicParticle * daughterparticle0 G4DynamicParticle * daughterparticle0
= new G4DynamicParticle( daughters[0], daughtermomentum[0]*direction0); = new G4DynamicParticle( G4MT_daughters[0], daughtermomentum[0]*direction0);
products->PushProducts(daughterparticle0); products->PushProducts(daughterparticle0);
@ -212,15 +211,15 @@ G4DecayProducts *MuDecayChannelWithSpin::DecayIt(G4double)
G4double beta = -1.0*daughtermomentum[0]/energy2; G4double beta = -1.0*daughtermomentum[0]/energy2;
G4double costhetan = 2.*G4UniformRand()-1.0; G4double costhetan = 2.*G4UniformRand()-1.0;
G4double sinthetan = std::sqrt((1.0-costhetan)*(1.0+costhetan)); G4double sinthetan = std::sqrt((1.0-costhetan)*(1.0+costhetan));
G4double phin = twopi*G4UniformRand()*rad; G4double phin = CLHEP::twopi*G4UniformRand()*CLHEP::rad;
G4double sinphin = std::sin(phin); G4double sinphin = std::sin(phin);
G4double cosphin = std::cos(phin); G4double cosphin = std::cos(phin);
G4ThreeVector direction1(sinthetan*cosphin,sinthetan*sinphin,costhetan); G4ThreeVector direction1(sinthetan*cosphin,sinthetan*sinphin,costhetan);
G4DynamicParticle * daughterparticle1 G4DynamicParticle * daughterparticle1
= new G4DynamicParticle( daughters[1], direction1*(vmass/2.)); = new G4DynamicParticle( G4MT_daughters[1], direction1*(vmass/2.));
G4DynamicParticle * daughterparticle2 G4DynamicParticle * daughterparticle2
= new G4DynamicParticle( daughters[2], direction1*(-1.0*vmass/2.)); = new G4DynamicParticle (G4MT_daughters[2], direction1*(-1.0*vmass/2.));
// boost to the muon rest frame // boost to the muon rest frame
G4LorentzVector p4; G4LorentzVector p4;

12
src/musrSteppingAction.cc
View File

@ -110,10 +110,10 @@ void musrSteppingAction::UserSteppingAction(const G4Step* aStep) {
sprintf(eMessage,"musrSteppingAction: Current number of steps for the track > %d ==> TRACK KILLED", sprintf(eMessage,"musrSteppingAction: Current number of steps for the track > %d ==> TRACK KILLED",
musrParameters::maximumNrOfStepsPerTrack); musrParameters::maximumNrOfStepsPerTrack);
musrErrorMessage::GetInstance()->musrError(WARNING,eMessage,true); musrErrorMessage::GetInstance()->musrError(WARNING,eMessage,true);
G4double x = postStepPosition.x()/mm; G4double x = postStepPosition.x()/CLHEP::mm;
G4double y = postStepPosition.y()/mm; G4double y = postStepPosition.y()/CLHEP::mm;
G4double z = postStepPosition.z()/mm; G4double z = postStepPosition.z()/CLHEP::mm;
G4double E=aTrack->GetVertexKineticEnergy()/MeV; G4double E=aTrack->GetVertexKineticEnergy()/CLHEP::MeV;
myRootOutput->htest1->Fill(x,y); myRootOutput->htest1->Fill(x,y);
myRootOutput->htest2->Fill(sqrt(x*x+y*y),z); myRootOutput->htest2->Fill(sqrt(x*x+y*y),z);
myRootOutput->htest3->Fill(E); myRootOutput->htest3->Fill(E);
@ -136,13 +136,13 @@ void musrSteppingAction::UserSteppingAction(const G4Step* aStep) {
// Temporary fix to avoid crashes caused by particles with unphysically high energies // Temporary fix to avoid crashes caused by particles with unphysically high energies
// (probably corrupted event?) // (probably corrupted event?)
// if ((aStep->GetPreStepPoint()->GetKineticEnergy()) > (1*GeV)) { // if ((aStep->GetPreStepPoint()->GetKineticEnergy()) > (1*GeV)) {
if ((preStepPoint->GetKineticEnergy()) > (1*GeV)) { if ((preStepPoint->GetKineticEnergy()) > (1*CLHEP::GeV)) {
musrErrorMessage::GetInstance()->musrError(SERIOUS, musrErrorMessage::GetInstance()->musrError(SERIOUS,
"musrSteppingAction: kinetic energy of a particle larger than 1GeV! STRANGE FOR muSR!",false); "musrSteppingAction: kinetic energy of a particle larger than 1GeV! STRANGE FOR muSR!",false);
G4RunManager* fRunManager = G4RunManager::GetRunManager(); G4RunManager* fRunManager = G4RunManager::GetRunManager();
G4cout<<" Event nr.:"<<fRunManager->GetCurrentEvent()->GetEventID() G4cout<<" Event nr.:"<<fRunManager->GetCurrentEvent()->GetEventID()
<<", the particle \""<< aTrack->GetDynamicParticle()->GetDefinition()->GetParticleName() <<", the particle \""<< aTrack->GetDynamicParticle()->GetDefinition()->GetParticleName()
<<"\" has energy of "<<(preStepPoint->GetKineticEnergy())/GeV<<" GeV."<<G4endl; <<"\" has energy of "<<(preStepPoint->GetKineticEnergy())/CLHEP::GeV<<" GeV."<<G4endl;
G4cout<<" Deleting the event!"<<G4endl; G4cout<<" Deleting the event!"<<G4endl;
G4cout.flush(); G4cout.flush();
myRootOutput->SetEventWeight(0); myRootOutput->SetEventWeight(0);