diff --git a/include/PrimaryGeneratorAction.hh b/include/PrimaryGeneratorAction.hh index 52baeaeae65ff15ac6ad64914b035eeead28860b..dbf111fee5c2b581288f5b7bb4532eb07acc0ca3 100644 --- a/include/PrimaryGeneratorAction.hh +++ b/include/PrimaryGeneratorAction.hh @@ -30,6 +30,7 @@ #ifndef PrimaryGeneratorAction_h #define PrimaryGeneratorAction_h 1 +#include "Randomize.hh" #include "G4VUserPrimaryGeneratorAction.hh" #include "globals.hh" #include <iostream> @@ -94,6 +95,8 @@ private: G4ParticleDefinition* fD0B; G4ParticleDefinition* fDsP; G4ParticleDefinition* fDsM; + + G4RandGeneral* rand_general; }; //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... diff --git a/run1.mac b/run1.mac index 75d9e09e19df596977b7348ebafc578d4e4066ad..889b6c7eab2da6d911531281ec0cd77263002309 100644 --- a/run1.mac +++ b/run1.mac @@ -10,6 +10,6 @@ #/process/list #/tracking/verbose 1 -/run/beamOn 10000 +/run/beamOn 10 #/run/particle/dumpList diff --git a/src/PrimaryGeneratorAction.cc b/src/PrimaryGeneratorAction.cc index a0870b5094198c8dcadf9c56dc02fa72792d0b51..327cf72ac7d1501b26c4b2407771c7e2d7a5847b 100644 --- a/src/PrimaryGeneratorAction.cc +++ b/src/PrimaryGeneratorAction.cc @@ -92,7 +92,8 @@ PrimaryGeneratorAction::PrimaryGeneratorAction() fD0(nullptr), fD0B(nullptr), fDsP(nullptr), - fDsM(nullptr) + fDsM(nullptr), + rand_general(nullptr) { G4AutoLock lock(&myMutex); @@ -139,7 +140,8 @@ PrimaryGeneratorAction::PrimaryGeneratorAction() fD0B = particleTable->FindParticle("anti_D0"); fDsP = particleTable->FindParticle("Ds+"); fDsM = particleTable->FindParticle("Ds-"); - + double prob[5] = {2241,217,56,5,2}; + rand_general = new G4RandGeneral(prob,5); if(!source) { //source = new ifstream("gen_ntuple2/input/numu_faser_1M.dump._001.txt"); source = new ifstream("../fnumin/electron.txt"); @@ -383,9 +385,20 @@ void PrimaryGeneratorAction::GeneratePrimaries(G4Event* anEvent) } //printf("Xin: %d %f %f %f %f\n",i1,tmp1,tmp2,tmp3,tmp4); + double rnd = rand_general->shoot(); + double_t tt = 0; + if(rnd<0.2) tt = G4RandFlat::shoot(0.01,0.2); + else if(rnd<0.2) tt = G4RandFlat::shoot(0.2,0.4); + else if(rnd<0.2) tt = G4RandFlat::shoot(0.4,0.6); + else if(rnd<0.2) tt = G4RandFlat::shoot(0.6,0.8); + else tt = G4RandFlat::shoot(0.8,1); + + + + fParticleGun->SetParticleEnergy(tmp1*GeV); - fParticleGun->SetParticleMomentumDirection(G4ThreeVector(0,0,1)); + fParticleGun->SetParticleMomentumDirection(G4ThreeVector(0,tt,1)); fParticleGun->SetParticlePosition(G4ThreeVector(0,0,-100.*cm)); fParticleGun->GeneratePrimaryVertex(anEvent); }