diff --git a/include/musrPrimaryGeneratorAction.hh b/include/musrPrimaryGeneratorAction.hh index 4965535..0906a86 100644 --- a/include/musrPrimaryGeneratorAction.hh +++ b/include/musrPrimaryGeneratorAction.hh @@ -65,6 +65,7 @@ class musrPrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction void SetTilt(G4ThreeVector v) {xangle0=v[0]; yangle0=v[1];} void SetSigmaTilt(G4ThreeVector v) {xangleSigma=v[0]; yangleSigma=v[1];} void SetPitch(G4double val) {pitch=val;} + void SetBeamDirection(G4ThreeVector vIniDir); void SetInitialMuonPolariz(G4ThreeVector vIniPol); void SetInitialPolarizFraction(G4double val) { if ((val>1.)||(val<-1.)) { @@ -111,6 +112,7 @@ class musrPrimaryGeneratorAction : public G4VUserPrimaryGeneratorAction G4double xangle0, yangle0, xangleSigma, yangleSigma, pitch; G4bool UnpolarisedMuonBeam, TransversalyUnpolarisedMuonBeam; G4double xPolarisIni, yPolarisIni, zPolarisIni; + G4double xDirection, yDirection, zDirection; G4double polarisFraction; G4double muonDecayTimeMin; G4double muonDecayTimeMax; diff --git a/include/musrPrimaryGeneratorMessenger.hh b/include/musrPrimaryGeneratorMessenger.hh index ae7af0e..92ac594 100644 --- a/include/musrPrimaryGeneratorMessenger.hh +++ b/include/musrPrimaryGeneratorMessenger.hh @@ -65,6 +65,7 @@ class musrPrimaryGeneratorMessenger: public G4UImessenger G4UIcmdWith3VectorAndUnit* setSigmaTiltAngleCmd; G4UIcmdWithADoubleAndUnit* setPitchCmd; G4UIcmdWith3Vector* setMuonPolarizCmd; + G4UIcmdWith3Vector* setDirectionCmd; G4UIcmdWithADouble* setMuonPolarizFractionCmd; G4UIcmdWith3VectorAndUnit* setMuonDecayTimeCmd; G4UIcmdWithAString* setTurtleCmd; diff --git a/src/musrPrimaryGeneratorAction.cc b/src/musrPrimaryGeneratorAction.cc index 6c066d9..6d79c57 100644 --- a/src/musrPrimaryGeneratorAction.cc +++ b/src/musrPrimaryGeneratorAction.cc @@ -271,8 +271,8 @@ void musrPrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) G4double px, py, pz; px = p*sin(xangle); py = p*sin(yangle); - pz = std::sqrt(p*p - px*px - py*py); - + // works for beam along +-z only for now + pz = zDirection * std::sqrt(p*p - px*px - py*py); // Assign spin G4double xpolaris=0, ypolaris=0, zpolaris=0; @@ -381,6 +381,16 @@ void musrPrimaryGeneratorAction::SetInitialMuonPolariz(G4ThreeVector vIniPol) } } +//=============================================================================== +void musrPrimaryGeneratorAction::SetBeamDirection(G4ThreeVector vIniDir) +{ + G4double magnitude=vIniDir.mag(); + xDirection=vIniDir(0)/magnitude; + yDirection=vIniDir(1)/magnitude; + zDirection=vIniDir(2)/magnitude; + G4cout<< "Initial Beam Direction set to ("<SetDefaultValue(0) ; setMuonPolarizCmd = new G4UIcmdWith3Vector("/gun/muonPolarizVector",this); - setMuonPolarizCmd->SetGuidance("Set initial mu polarisation as a vector (without unit)"); + setMuonPolarizCmd->SetGuidance("Set initial mu polarisation as a vector (without units)"); setMuonPolarizCmd->SetGuidance(" The vector does not have to be normalised to 1"); setMuonPolarizCmd->SetParameterName("mes_polarisX","mes_polarisY","mes_polarisZ",true,true); + setDirectionCmd = new G4UIcmdWith3Vector("/gun/direction",this); + setDirectionCmd->SetGuidance("Set initial mu beam direction as a vector (without units)"); + setDirectionCmd->SetGuidance(" The vector does not have to be normalised to 1"); + setDirectionCmd->SetParameterName("DirectionX","DirectionY","DirectionZ",true,true); + setMuonPolarizFractionCmd = new G4UIcmdWithADouble("/gun/muonPolarizFraction",this); setMuonPolarizFractionCmd->SetGuidance(" Set the fraction of the muon polarisation (in the range of -1 to 1),"); setMuonPolarizFractionCmd->SetGuidance(" where fraction = (N_up_spin - N_down_spin) / (N_up_spin + N_down_spin)"); @@ -196,6 +201,7 @@ musrPrimaryGeneratorMessenger::~musrPrimaryGeneratorMessenger() delete setSigmaTiltAngleCmd; delete setPitchCmd; delete setMuonPolarizCmd; + delete setDirectionCmd; delete setMuonPolarizFractionCmd; delete setMuonDecayTimeCmd; delete setTurtleCmd; @@ -245,6 +251,8 @@ void musrPrimaryGeneratorMessenger::SetNewValue(G4UIcommand * command,G4String n { musrAction->SetPitch(setPitchCmd->GetNewDoubleValue(newValue));} if( command == setMuonPolarizCmd) { musrAction->SetInitialMuonPolariz(setMuonPolarizCmd->GetNew3VectorValue(newValue));} + if( command == setDirectionCmd) + { musrAction->SetBeamDirection(setDirectionCmd->GetNew3VectorValue(newValue));} if( command == setMuonPolarizFractionCmd) { musrAction->SetInitialPolarizFraction(setMuonPolarizFractionCmd->GetNewDoubleValue(newValue));} if( command == setMuonDecayTimeCmd)