Commit 85673e78 authored by Elidiano Tronchin's avatar Elidiano Tronchin
Browse files

* Chqnge to generate muon

parent d2ba76e7
......@@ -240,6 +240,8 @@ string ParticleStack::ChangeNameToG4Standard(string OriginalName){
else if (FinalName=="deuteron") FinalName="deuteron";
else if (FinalName=="triton") FinalName="triton";
else if (FinalName=="alpha") FinalName="alpha";
else if (FinalName=="mu+") FinalName="mu+";
else if (FinalName=="mu-") FinalName="mu-";
else if (FinalName=="n") FinalName="neutron";
else if (FinalName=="neutron") FinalName="neutron";
return(FinalName);
......
......@@ -368,7 +368,6 @@ void Dali::ReadSensitive(const G4Event* ){
unsigned int size = Scorer->GetMult();
// cout << "size " << size << endl;
for(unsigned int i = 0 ; i < size ; i++){
vector<unsigned int> level = Scorer->GetLevel(i);
double Energy = RandGauss::shoot(Scorer->GetEnergy(i),Dali_NS::ResoEnergy);
......
......@@ -174,8 +174,8 @@ void EventGeneratorCosmic::GenerateEvent(G4Event*){
/* //Putting a cylinder
if(randomize>0){ //top
momentum_x = cos(angle)*CosmicAngle;
momentum_z = sin(angle)*CosmicAngle;
momentum_x = cos(angle)*sin(CosmicAngle);
momentum_z = sin(angle)*sin(CosmicAngle);
x0 = cos(angle2*2)*R*(randomize*2);
z0 = sin(angle2*2)*R*(randomize*2);
......@@ -186,17 +186,19 @@ void EventGeneratorCosmic::GenerateEvent(G4Event*){
z0 = sin(angle)*R;
par.m_y0 = (.5-RandFlat::shoot())*H; ///!!!!!!!
momentum_x = -cos(angle+angle2)*CosmicAngle;
momentum_z = -sin(angle+angle2)*CosmicAngle;
momentum_x = -cos(angle+angle2)*sin(CosmicAngle);
momentum_z = -sin(angle+angle2)*sin(CosmicAngle);
}
*/
// Begin Constrain to pass in a circle with radius 3R
momentum_x = cos(angle)*CosmicAngle;
momentum_z = sin(angle)*CosmicAngle;
momentum_y = cos(CosmicAngle);
momentum_x = cos(angle)*sin(CosmicAngle);
momentum_z = sin(angle)*sin(CosmicAngle);
x0 = cos(angle2*2)*R*(0.5-randomize)*3-momentum_x*( H/2 / momentum_y);// *( H/2 / momentum_y) this is to have the origin always with par.m_y0 = H/2;
z0 = sin(angle2*2)*R*(0.5-randomize)*3-momentum_z*( H/2 / momentum_y);
par.m_y0 = H/2;
par.m_y0 = H/2; // momentum_y*( H/2 / momentum_y);
// End Constrain to pass in a circle with radius 3R
......
......@@ -91,6 +91,8 @@ void EventGeneratorIsotropic::ReadConfiguration(NPL::InputParser parser){
else if(particleName[j]=="3He" || particleName[j]=="He3") { it->m_particleName.push_back("3He") ; }
else if(particleName[j]=="alpha") { it->m_particleName.push_back("4He") ; }
else if(particleName[j]=="gamma") { it->m_particleName.push_back("gamma") ;}
else if(particleName[j]=="mu+") { it->m_particleName.push_back("mu+") ;}
else if(particleName[j]=="mu-") { it->m_particleName.push_back("mu-") ;}
else if(particleName[j]=="neutron") {it->m_particleName.push_back("neutron") ;}
else it->m_particleName.push_back(particleName[j]);
}
......@@ -128,7 +130,7 @@ void EventGeneratorIsotropic::GenerateEvent(G4Event*){
par.m_particle=NULL;
if(par.m_particle==NULL){
if(par.m_particleName[p]=="gamma" || par.m_particleName[p]=="neutron" || par.m_particleName[p]=="opticalphoton" || par.m_particleName[p]=="mu+"){
if(par.m_particleName[p]=="gamma" || par.m_particleName[p]=="neutron" || par.m_particleName[p]=="opticalphoton" || par.m_particleName[p]=="mu+" || par.m_particleName[p]=="mu-"){
par.m_particle = G4ParticleTable::GetParticleTable()->FindParticle(par.m_particleName[p].c_str());
}
else{
......
......@@ -118,6 +118,7 @@ void EventGeneratorMultipleParticle::ReadConfiguration(NPL::InputParser parser){
else if(sParticle=="3He" || sParticle=="He3"){vParticle.push_back("3He");}
else if(sParticle=="alpha"){vParticle.push_back("4He");}
else if(sParticle=="gamma") { vParticle.push_back("gamma") ;}
else if(sParticle=="mu+") { vParticle.push_back("mu+") ;}
else if(sParticle=="neutron") {vParticle.push_back("neutron") ;}
else vParticle.push_back(sParticle);
......@@ -154,7 +155,7 @@ void EventGeneratorMultipleParticle::GenerateEvent(G4Event* evt){
for(int i=0; i<m_Multiplicity[evtID]; i++){
m_particle=NULL;
if(m_particle==NULL){
if(m_particleName[evtID][i]=="gamma" || m_particleName[evtID][i]=="neutron" || m_particleName[evtID][i]=="opticalphoton"){
if(m_particleName[evtID][i]=="gamma" || m_particleName[evtID][i]=="neutron" || m_particleName[evtID][i]=="opticalphoton" || m_particleName[evtID][i]=="mu+"){
m_particle = G4ParticleTable::GetParticleTable()->FindParticle(m_particleName[evtID][i].c_str());
}
else{
......
......@@ -43,6 +43,7 @@
#include "G4IonFluctuations.hh"
#include "G4IonParametrisedLossModel.hh"
#include "G4UniversalFluctuation.hh"
#include "G4UImanager.hh"
#include "G4PAIModel.hh"
#include "G4PAIPhotModel.hh"
......@@ -50,6 +51,11 @@
/////////////////////////////////////////////////////////////////////////////
PhysicsList::PhysicsList() : G4VUserPhysicsList(){
// UI commands to activate one or more High Energy Processes
// G4UImanager* UI = G4UImanager::GetUIpointer();
// UI->ApplyCommand("/physics_list/em/MuonNuclear true");
m_EmList = "Option4";
defaultCutValue = 1*mm;//0.2*mm;
opticalPhysicsList = NULL;
......
......@@ -6,6 +6,8 @@
#include "G4PhysListFactory.hh"
// UI
#include "G4UImanager.hh"
#include "QBBC.hh"
#include "G4UIterminal.hh"
#include "G4UItcsh.hh"
......@@ -72,7 +74,14 @@ int main(int argc, char** argv){
runManager->SetUserInitialization(detector);
PhysicsList* physicsList = new PhysicsList();
//G4VModularPhysicsList* physicsList = new QBBC;
physicsList->SetVerboseLevel(0);
runManager->SetUserInitialization(physicsList);
PrimaryGeneratorAction* primary = new PrimaryGeneratorAction(detector);
// Initialize Geant4 kernel
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment