Difference between revisions of "DV Creating LUND Files"
Line 205: | Line 205: | ||
=Writing LUND files= | =Writing LUND files= | ||
+ | <pre> | ||
+ | void MollerLUND() | ||
+ | { | ||
+ | struct evt_t | ||
+ | { | ||
+ | Long_t event; | ||
+ | Float_t IntKE, IntMom[3],IntPos[3],FnlKE,FnlMom[3],FnlPos[3],MolKE,MolMom[3],MolPos[3]; | ||
+ | }; | ||
+ | ifstream in; | ||
+ | in.open("MollerScattering_NH3_4e8incident.dat"); | ||
+ | evt_t evt; | ||
+ | Int_t nlines=0; | ||
+ | TFile *f = new TFile("MollerScattering_NH3_4e8incident.root","RECREATE"); | ||
+ | TTree *tree0 = new TTree("MollerScattering_NH3_4e8","Moller data from ascii file"); | ||
+ | FILE *f0; | ||
+ | f0=fopen("MollerScattering_NH3_4e8.LUND","w"); | ||
+ | |||
+ | tree0->Branch("evt",&evt.event,"event/I:IntKE/F:IntPx:IntPy:IntPz:IntPosx:IntPosy:IntPosz:FnlKE:FnlPx:FnlPy:FnlPz:FnlPosx:FnlPosy:FnlPosz: | ||
+ | MolKE:MolPx:MolPy:MolPz:MolPosx:MolPosy:MolPosz"); | ||
+ | |||
+ | while(in.good()) | ||
+ | { | ||
+ | evt.event=nlines; | ||
+ | in >> evt.IntKE >> evt.IntMom[0] >> evt.IntMom[1] >> evt.IntMom[2] >> evt.IntPos[0] | ||
+ | >> evt.IntPos[1] >> evt.IntPos[2] >> evt.FnlKE >> evt.FnlMom[0] >> evt.FnlMom[1] >> evt.FnlMom[2] | ||
+ | >> evt.FnlPos[0] >> evt.FnlPos[1] >> evt.FnlPos[2] >> evt.MolKE >> evt.MolMom[0] >> evt.MolMom[1] | ||
+ | >> evt.MolMom[2] >> evt.MolPos[0] >> evt.MolPos[1] >> evt.MolPos[2]; double Nu; | ||
+ | double x_bj; // Bjorken-x | ||
+ | double y; | ||
+ | double W ; // missing mass | ||
+ | |||
+ | // Momentum and Energy are measured in MeV! | ||
+ | |||
+ | evt.FnlKE=sqrt(evt.FnlMom[0]*evt.FnlMom[0]+evt.FnlMom[1]*evt.FnlMom[1]+evt.FnlMom[2]*evt.FnlMom[2]+0.511*0.511); | ||
+ | evt.MolKE=sqrt(evt.MolMom[0]*evt.MolMom[0]+evt.MolMom[1]*evt.MolMom[1]+evt.MolMom[2]*evt.MolMom[2]+0.511*0.511); | ||
+ | |||
+ | |||
+ | Nu=evt.IntKE-evt.FnlKE; | ||
+ | Qsqrd=4*evt.IntKE*evt.FnlKE*(1-evt.FnlMom[2]/evt.FnlKE)/2; /* should be final momentum and not final energy*/ | ||
+ | |||
+ | x_bj=Qsqrd/(2*0.938*Nu); | ||
+ | y=Nu/evt.IntKE; | ||
+ | W=0.938*0.938+2*0.938*Nu-Qsqrd; | ||
+ | if(W>0) | ||
+ | W=sqrt(W); | ||
+ | |||
+ | // | ||
+ | // MeV MUST be converted to GeV for LUND format! | ||
+ | |||
+ | double px,py,pz; | ||
+ | double Px,Py,Pz; | ||
+ | double KE,ke; | ||
+ | |||
+ | Px=evt.FnlMom[0]/10; | ||
+ | Py=evt.FnlMom[1]/10; | ||
+ | Pz=evt.FnlMom[2]/10; | ||
+ | px=evt.MolMom[0]/100; | ||
+ | py=evt.MolMom[1]/100; | ||
+ | pz=evt.MolMom[2]/100; | ||
+ | |||
+ | KE=evt.FnlKE/10; | ||
+ | ke=evt.MolKE/100; | ||
+ | |||
+ | |||
+ | fprintf(f0,"%d\t%g\t%g\t%g\t%g\t%g\t%g\t%g\t%g\t%g\n",2,1.,1.,1.,1.,x_bj,y,W,Qsqrd,Nu); | ||
+ | fprintf(f0, "%d\t%g\t%g\t%g\t%g\t%g\t%f\t%f\t%f\t%f\t%g\t%g\t%g\t%g\n", | ||
+ | 1,-1.,1.,11.,0.,0.,Px,Py,Pz,KE,0.000511, 0.,0.,0.); | ||
+ | fprintf(f0, "%d\t%g\t%g\t%g\t%g\t%g\t%f\t%f\t%f\t%f\t%g\t%g\t%g\t%g\n", | ||
+ | 2,-1.,1.,11.,0.,0.,px,py,pz,ke,0.000511, 0.,0.,0.); | ||
+ | |||
+ | |||
+ | nlines++; | ||
+ | tree0->Fill(); | ||
+ | } | ||
+ | printf("NumberLines=%d\n",nlines); | ||
+ | tree0->Print(); | ||
+ | tree0->Write(); | ||
+ | delete tree0; | ||
+ | |||
+ | } | ||
+ | </pre> | ||
=Correct for vertex location= | =Correct for vertex location= |
Revision as of 20:48, 29 March 2016
The LUND format
The LUND file format is broken into two parts. The first part of the format is the header, which basically tells how many particles follow this line. For Moller scattering, this number should always be two electrons. The second component of the LUND format contains the kinematic variable for the scattered electron and the Moller electron.
The LUND format has extra variables that are not utilized within GEMC. Only the BOLD variables are necessary within GEMC simulations.
The Header
Column | Quantity |
---|---|
Where
1:Number of Particles
This line tells how many particles follow the header line. For Moller Scattering, this number should always be 2.
2:Number of Target Nucleons
For this simulation, only an electron-electron collision is considered. This quantity is always set to 1 for the one stationary electron we consider the incident electron scattering from.
3:Number of Target Protons
For Moller Scattering, there are no target protons. This number is set to 1, but does not have any effect within the GEMC simulations.
4:Target Polarization
This represents the polarization of the target material, either positive or negative 1. This value is always set to 1 and has no effect within the GEMC simulations.
5:Beam Polarization
This represents the polarization of the electron beam, either negative or positive 1. This value is always set to positive 1.
6:x
This represents the Bjorken x scaling variable.
Where M is the rest mass-energy of the proton
7:y
8:W
The invariant mass of the final hadronic system
9:
This represents the squared 4-momentum-transfer vector q of the exchanged virtual photon.
Where q is the momentum transfer betwwn the incident electron and target via the virtual photon.
10:
This represents the energy loss between scattering electrons.
This can be written in the Lab frame as:
where Ei and Ef are the initial and final electron energies.
Particles
Column | Quantity |
---|---|
Where
1:Index
From the header line 1, which states how many particles follow, this counts in increments of one up to the required number.
2:Charge
-1, 0, 1
3:Type
1 denotes that this is an active particle that can interact within the simulation.
4:Particle ID
This quantity follows the Particle Data Group Monte Carlo numbering scheme for representing what type of particle is simulated.
File:MontecarlorPDGparticleID.pdf
5:Parent Index
This signifies if the particle can be considered a primary (=0) or secondary particle. For this simulation, GEMC does not utilize this variable and all values have been set to 0
6:Daughter Index
This signifies if the particle has created secondary particles. For this simulation, GEMC does not utilize this variable and all values have been set to 0.
7:Momentum x
This is the value of the momentum vector in the x direction. Values must be given in GeV.
8:Momentum y
This is the value of the momentum vector in the y direction. Values must be given in GeV.
9:Momentum z
This is the value of the momentum vector in the z direction. Values must be given in GeV.
10:E
E is the total relativistic energy,
For the Ultrarelativistic limit,
Since
11:Mass
This represents the rest mass-energy of the particle being simulated. For Moller scattering this value has been set to .000511 GeV.
12:Vertex x
This represents the x position within the eg12 detector at time t=0.
13:Vertex y
This represents the y position within the eg12 detector at time t=0.
14:Vertex z
This represents the z position within the eg12 detector at time t=0.
Writing LUND files
void MollerLUND() { struct evt_t { Long_t event; Float_t IntKE, IntMom[3],IntPos[3],FnlKE,FnlMom[3],FnlPos[3],MolKE,MolMom[3],MolPos[3]; }; ifstream in; in.open("MollerScattering_NH3_4e8incident.dat"); evt_t evt; Int_t nlines=0; TFile *f = new TFile("MollerScattering_NH3_4e8incident.root","RECREATE"); TTree *tree0 = new TTree("MollerScattering_NH3_4e8","Moller data from ascii file"); FILE *f0; f0=fopen("MollerScattering_NH3_4e8.LUND","w"); tree0->Branch("evt",&evt.event,"event/I:IntKE/F:IntPx:IntPy:IntPz:IntPosx:IntPosy:IntPosz:FnlKE:FnlPx:FnlPy:FnlPz:FnlPosx:FnlPosy:FnlPosz: MolKE:MolPx:MolPy:MolPz:MolPosx:MolPosy:MolPosz"); while(in.good()) { evt.event=nlines; in >> evt.IntKE >> evt.IntMom[0] >> evt.IntMom[1] >> evt.IntMom[2] >> evt.IntPos[0] >> evt.IntPos[1] >> evt.IntPos[2] >> evt.FnlKE >> evt.FnlMom[0] >> evt.FnlMom[1] >> evt.FnlMom[2] >> evt.FnlPos[0] >> evt.FnlPos[1] >> evt.FnlPos[2] >> evt.MolKE >> evt.MolMom[0] >> evt.MolMom[1] >> evt.MolMom[2] >> evt.MolPos[0] >> evt.MolPos[1] >> evt.MolPos[2]; double Nu; double x_bj; // Bjorken-x double y; double W ; // missing mass // Momentum and Energy are measured in MeV! evt.FnlKE=sqrt(evt.FnlMom[0]*evt.FnlMom[0]+evt.FnlMom[1]*evt.FnlMom[1]+evt.FnlMom[2]*evt.FnlMom[2]+0.511*0.511); evt.MolKE=sqrt(evt.MolMom[0]*evt.MolMom[0]+evt.MolMom[1]*evt.MolMom[1]+evt.MolMom[2]*evt.MolMom[2]+0.511*0.511); Nu=evt.IntKE-evt.FnlKE; Qsqrd=4*evt.IntKE*evt.FnlKE*(1-evt.FnlMom[2]/evt.FnlKE)/2; /* should be final momentum and not final energy*/ x_bj=Qsqrd/(2*0.938*Nu); y=Nu/evt.IntKE; W=0.938*0.938+2*0.938*Nu-Qsqrd; if(W>0) W=sqrt(W); // // MeV MUST be converted to GeV for LUND format! double px,py,pz; double Px,Py,Pz; double KE,ke; Px=evt.FnlMom[0]/10; Py=evt.FnlMom[1]/10; Pz=evt.FnlMom[2]/10; px=evt.MolMom[0]/100; py=evt.MolMom[1]/100; pz=evt.MolMom[2]/100; KE=evt.FnlKE/10; ke=evt.MolKE/100; fprintf(f0,"%d\t%g\t%g\t%g\t%g\t%g\t%g\t%g\t%g\t%g\n",2,1.,1.,1.,1.,x_bj,y,W,Qsqrd,Nu); fprintf(f0, "%d\t%g\t%g\t%g\t%g\t%g\t%f\t%f\t%f\t%f\t%g\t%g\t%g\t%g\n", 1,-1.,1.,11.,0.,0.,Px,Py,Pz,KE,0.000511, 0.,0.,0.); fprintf(f0, "%d\t%g\t%g\t%g\t%g\t%g\t%f\t%f\t%f\t%f\t%g\t%g\t%g\t%g\n", 2,-1.,1.,11.,0.,0.,px,py,pz,ke,0.000511, 0.,0.,0.); nlines++; tree0->Fill(); } printf("NumberLines=%d\n",nlines); tree0->Print(); tree0->Write(); delete tree0; }