Vanwdani: /* Writing LUND files */

2016-04-05T20:42:34Z

Writing LUND files

Vanwdani: /* References */

2016-04-05T19:57:41Z

References

Vanwdani: /* Links */

2016-04-05T19:57:31Z

Links

Vanwdani: /* Writing LUND files */

2016-04-05T19:57:10Z

Writing LUND files

Vanwdani: Created page with "=Writing LUND files= Reading the GEANT4 simulation output in,

 void MollerLUND() {         struct evt_t                 {                 Long_t event;                 Float…"

2016-04-05T19:55:19Z

Created page with "=Writing LUND files= Reading the GEANT4 simulation output in, <pre> void MollerLUND() { struct evt_t { Long_t event; Float…"

New page

=Writing LUND files=
Reading the GEANT4 simulation output in,
<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_4e7incident.root","RECREATE");
TTree *tree0 = new TTree("MollerScattering_NH3_4e7","Moller data from ascii file");
FILE *f0;
f0=fopen("MollerScattering_NH3_4e7.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>

The file is written as seen below.

[[File:LUND_view.png]]

=Corrections=
The prefered coordinates for the reconstruction within GEMC would have the vertex location at (0, 0, 0). The Moller event data from GEANT4 contains the location within that simulation where the scattering occured. These variables will need to be renormalized within the LUND file to zero. Only the 14th column within the 2nd and 3rd rows out of every 3 rows needs to be changed. Using a script:

<pre>awk '{if('NR%3==2'||'NR%3==0') {$14="0"}; print $0}' MollerScattering_NH3_4e7.LUND>MollerScattering.LUND</pre>

The LUND format has a limit of 75000 lines, so to account for this by creating 21 new files with the suffix of aa-au are created:

<pre>split -l 75000 MollerScattering.LUND MollerScattering_NH3_4e7_ </pre>

We will use 5 of the 21 sets of 75000 entries in GEMC trials involving different solenoid field strengths. We can use the simple command:

<pre> mv MollerScattering_NH3_4e7_aa MollerScattering_NH3_4e7_aa.LUND</pre>

Repeating this for the other 4 files.

=References=
[1] M. Ungaro, ''et al.,'' Jefferson Lab:The LUND format, https://gemc.jlab.org/gemc/html/documentation/generator/lund.html (2016)

[2] F. Halzen, and A. Martin, ''Quarks and Leptons:An Introductory Course in Modern Particle Physics,'' John Wiley & Sons (1984)

[3] W. Tung, ''Bjorken scaling,'' http://www.scholarpedia.org/article/Bjorken_scaling (2009)

[4] L. Garren, ''et al.,'' MONTE CARLO PARTICLE NUMBERING SCHEME, http://pdg.lbl.gov/2007/reviews/montecarlorpp.pdf (2007)

----

← Older revision		Revision as of 19:57, 5 April 2016
Line 116:		Line 116:

	----		----
		+	=Links=
		+
		+	[[DV_RunGroupC_Moller#Preparing_LUND_for_Moller_Distribution\|Back]]

← Older revision		Revision as of 19:57, 5 April 2016
Line 89:		Line 89:

	----		----
−
−	~~=Links=~~
−
−	~~[[DV_RunGroupC_Moller#Preparing_LUND_for_Moller_Distribution\|Back]]~~

	=Corrections=		=Corrections=

← Older revision		Revision as of 19:57, 5 April 2016
Line 87:		Line 87:

	[[File:LUND_view.png]]		[[File:LUND_view.png]]
		+
		+	----
		+
		+	=Links=
		+
		+	[[DV_RunGroupC_Moller#Preparing_LUND_for_Moller_Distribution\|Back]]

	=Corrections=		=Corrections=

DV Running LUND for Moller Distribution - Revision history

Vanwdani: /* Writing LUND files */

Vanwdani: /* References */

Vanwdani: /* Links */

Vanwdani: /* Writing LUND files */

Vanwdani: Created page with "=Writing LUND files= Reading the GEANT4 simulation output in,
void MollerLUND() { struct evt_t { Long_t event; Float…"

@@ Line 54: / Line 54: @@
                  double KE,ke;
-                 Px=evt.FnlMom[0]/10;
+                 Px=evt.FnlMom[0]/1000;
-                 Py=evt.FnlMom[1]/10;
+                 Py=evt.FnlMom[1]/1000;
-                 Pz=evt.FnlMom[2]/10;
+                 Pz=evt.FnlMom[2]/1000;
-                 px=evt.MolMom[0]/100;
+                 px=evt.MolMom[0]/1000;
-                 py=evt.MolMom[1]/100;
+                 py=evt.MolMom[1]/1000;
-                 pz=evt.MolMom[2]/100;
+                 pz=evt.MolMom[2]/1000;
-                 KE=evt.FnlKE/10;
+                 KE=evt.FnlKE/1000;
-                 ke=evt.MolKE/100;
+                 ke=evt.MolKE/1000;