BOSS 7.1.1
BESIII Offline Software System
Loading...
Searching...
No Matches
Macros | Functions
phokhara.cc File Reference
#include "Phokhara/PhokharaDef.h"
#include <fstream>
#include <iostream>
#include <iomanip>
#include "Phokhara/ranlxd.h"

Go to the source code of this file.

Macros

#define RLXDRESETF(SEED)
 
#define INPUT(NGES, NM, OUTFILE)
 
#define INITHISTO()
 
#define ENDHISTO()
 
#define WRITEEVENT()
 
#define RANLXDF(AR, VAL)
 
#define GEN_1PH(I, QQMIN, QQMAX, COS1MIN, COS1MAX, COS3MIN, COS3MAX)
 
#define GEN_2PH(I, QQMIN, COS1MIN, COS1MAX, COS2MIN, COS2MAX, COS3MIN, COS3MAX)
 

Functions

 PROTOCCALLSFSUB1 (RLXDRESETF, rlxdresetf, INTV)
 
 PROTOCCALLSFSUB3 (INPUT, input, PLONG, PINT, PSTRING)
 
 PROTOCCALLSFSUB0 (INITHISTO, inithisto)
 
 PROTOCCALLSFSUB0 (ENDHISTO, endhisto)
 
 PROTOCCALLSFSUB0 (WRITEEVENT, writeevent)
 
 PROTOCCALLSFSUB2 (RANLXDF, ranlxdf, DOUBLEV, INT)
 
 PROTOCCALLSFSUB7 (GEN_1PH, gen_1ph_, INT, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE)
 
 PROTOCCALLSFSUB8 (GEN_2PH, gen_2ph_, INT, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE)
 
int main ()
 

Macro Definition Documentation

◆ ENDHISTO

#define ENDHISTO ( )
Value:
CCALLSFSUB0(ENDHISTO,endhisto)
CCALLSFSUB0
#define CCALLSFSUB0(UN, LN)
Definition cfortran.h:1155
ENDHISTO
#define ENDHISTO()
Definition phokhara.cc:21

Definition at line 21 of file phokhara.cc.

Referenced by main().

◆ GEN_1PH

#define GEN_1PH ( I,
QQMIN,
QQMAX,
COS1MIN,
COS1MAX,
COS3MIN,
COS3MAX )
Value:
CCALLSFSUB7(GEN_1PH,gen_1ph_,INT, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE,I,QQMIN,QQMAX,COS1MIN,COS1MAX,COS3MIN,COS3MAX)
I
const DifComplex I
CCALLSFSUB7
#define CCALLSFSUB7(UN, LN, T1, T2, T3, T4, T5, T6, T7, A1, A2, A3, A4, A5, A6, A7)
Definition cfortran.h:1118
GEN_1PH
#define GEN_1PH(I, QQMIN, QQMAX, COS1MIN, COS1MAX, COS3MIN, COS3MAX)
Definition phokhara.cc:30

Definition at line 30 of file phokhara.cc.

Referenced by main().

◆ GEN_2PH

#define GEN_2PH ( I,
QQMIN,
COS1MIN,
COS1MAX,
COS2MIN,
COS2MAX,
COS3MIN,
COS3MAX )
Value:
CCALLSFSUB8(GEN_2PH,gen_2ph_,INT, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE, PDOUBLE,I,QQMIN,COS1MIN,COS1MAX,COS2MIN,COS2MAX,COS3MIN,COS3MAX)
CCALLSFSUB8
#define CCALLSFSUB8(UN, LN, T1, T2, T3, T4, T5, T6, T7, T8, A1, A2, A3, A4, A5, A6, A7, A8)
Definition cfortran.h:1120
GEN_2PH
#define GEN_2PH(I, QQMIN, COS1MIN, COS1MAX, COS2MIN, COS2MAX, COS3MIN, COS3MAX)
Definition phokhara.cc:33

Definition at line 33 of file phokhara.cc.

Referenced by main().

◆ INITHISTO

#define INITHISTO ( )
Value:
CCALLSFSUB0(INITHISTO,inithisto)
INITHISTO
#define INITHISTO()
Definition phokhara.cc:18

Definition at line 18 of file phokhara.cc.

Referenced by main().

◆ INPUT

#define INPUT ( NGES,
NM,
OUTFILE )
Value:
CCALLSFSUB3(INPUT,input,PLONG,PINT,PSTRING,NGES,NM,OUTFILE)
CCALLSFSUB3
#define CCALLSFSUB3(UN, LN, T1, T2, T3, A1, A2, A3)
Definition cfortran.h:1110
INPUT
#define INPUT(NGES, NM, OUTFILE)
Definition phokhara.cc:15

Definition at line 15 of file phokhara.cc.

Referenced by main().

◆ RANLXDF

#define RANLXDF ( AR,
VAL )
Value:
CCALLSFSUB2(RANLXDF,ranlxdf,DOUBLEV, INT, AR, VAL)
CCALLSFSUB2
#define CCALLSFSUB2(UN, LN, T1, T2, A1, A2)
Definition cfortran.h:1108
RANLXDF
#define RANLXDF(AR, VAL)
Definition phokhara.cc:27

Definition at line 27 of file phokhara.cc.

Referenced by main().

◆ RLXDRESETF

#define RLXDRESETF ( SEED)
Value:
CCALLSFSUB1(RLXDRESETF,rlxdresetf,INTV, SEED)
CCALLSFSUB1
#define CCALLSFSUB1(UN, LN, T1, A1)
Definition cfortran.h:1106
RLXDRESETF
#define RLXDRESETF(SEED)
Definition phokhara.cc:12

Definition at line 12 of file phokhara.cc.

Referenced by main().

◆ WRITEEVENT

#define WRITEEVENT ( )
Value:
CCALLSFSUB0(WRITEEVENT,writeevent)
WRITEEVENT
#define WRITEEVENT()
Definition phokhara.cc:24

Definition at line 24 of file phokhara.cc.

Function Documentation

◆ main()

int main ( )

Definition at line 35 of file phokhara.cc.

35 {
36 double qqmin,qqmax;
37 double cos1min,cos1max,cos2min,cos2max,cos3min,cos3max;
38 double dsigm1,dsigm2,sigma1,sigma2,sigma,dsigm,Ar[14],Ar_r[14];
39 int nm,i,s_seed[105];
40 long int nges,k,j;
41 char outfile[20];
42
43 // --- reads the seed ------
44 ifstream seeds("seed.dat");
45 if( seeds.is_open() )
46 {
47 int ii=0;
48 while( ! seeds.eof() )
49 seeds >> s_seed[ii++];
50 }
51 else
52 cerr << "Cannot open file seed.dat for reading" << endl;
53
54 RLXDRESETF(s_seed);
55 //rlxd_reset(s_seed);
56
57 // --- input parameters ----------------------------
58
59 INPUT(nges, nm, outfile);
60
61 // --- open output file for generated momenta ------
62 // if(iprint.ne.0) open (10,file=outfile,type="new")
63
64 // --- print run data ------------------------------
65 cout << "----------------------------------------------------" << FLAGS.pion<< endl;
66 if (FLAGS.pion == 0)
67 cout << " PHOKHARA 6.0 : e^+ e^- -> mu^+ mu^- gamma" << endl;
68 else if (FLAGS.pion == 1)
69 cout << " PHOKHARA 6.0: e^+ e^- -> pi^+ pi^- gamma" << endl;
70 else if (FLAGS.pion == 2)
71 cout << " PHOKHARA 6.0: e^+ e^- -> pi^+ pi^- 2pi^0 gamma" << endl;
72 else if (FLAGS.pion == 3)
73 cout << " PHOKHARA 6.0: e^+ e^- -> 2pi^+ 2pi^- gamma" << endl;
74 else if (FLAGS.pion == 4)
75 cout << " PHOKHARA 6.0: e^+ e^- -> p pbar gamma" << endl;
76 else if (FLAGS.pion == 5)
77 cout << " PHOKHARA 6.0: e^+ e^- -> n nbar gamma" << endl;
78 else if (FLAGS.pion == 6)
79 cout << " PHOKHARA 6.0: e^+ e^- -> K^+ K^- gamma" << endl;
80 else if (FLAGS.pion == 7)
81 cout << " PHOKHARA 6.0: e^+ e^- -> K_0 K_0bar gamma" << endl;
82 else if (FLAGS.pion == 8)
83 cout << " PHOKHARA 6.0: e^+ e^- -> pi^+ pi^- pi^0 gamma" << endl;
84 else if (FLAGS.pion == 9) {
85 cout << "PHOKHARA 6.0 : e^+ e^- ->" << endl;
86 cout << " Lambda (-> pi^- p) Lambda bar (-> pi^+ pbar) gamma" << endl;
87 } else
88 cout << " PHOKHARA 6.0: not yet implemented" << endl;
89
90 // --------------------------------
91 cout << "----------------------------------------------------" << endl;
92 printf(" %s %f %s\n","cms total energy = ",sqrt(CTES.Sp)," GeV ");
93 if (FLAGS.tagged == 0) {
94 if((CUTS.gmin/2.0/CTES.ebeam) < 0.0098){
95 cerr << " minimal missing energy set to small" << endl;
96 return 0;
97 }
98 printf(" %s %f %s\n","minimal tagged photon energy = ",CUTS.gmin," GeV ");
99 printf(" %s %f,%f\n","angular cuts on tagged photon = ",CUTS.phot1cut,CUTS.phot2cut);
100 }
101 else {
102 if((CUTS.gmin/2.0/CTES.ebeam) < 0.0098){
103 cerr << " minimal missing energy set to small" << endl;
104 return 0;
105 }
106 printf(" %s %f %s\n","minimal missing energy = ",CUTS.gmin," GeV ");
107 printf(" %s %f,%f\n","angular cuts on missing momentum = ",CUTS.phot1cut, CUTS.phot2cut);
108 }
109
110
111 // --------------------------------
112 if (FLAGS.pion == 0)
113 printf(" %s %f,%f\n","angular cuts on muons = ",CUTS.pi1cut,CUTS.pi2cut);
114 else if (FLAGS.pion == 4)
115 printf(" %s %f,%f\n","angular cuts on protons = ",CUTS.pi1cut,CUTS.pi2cut);
116 else if (FLAGS.pion == 5)
117 printf(" %s %f,%f\n","angular cuts on neutrons = ", CUTS.pi1cut,CUTS.pi2cut);
118 else if ((FLAGS.pion == 6)||(FLAGS.pion == 7))
119 printf(" %s %f,%f\n","angular cuts on kaons = ", CUTS.pi1cut,CUTS.pi2cut);
120 else if (FLAGS.pion == 9)
121 printf(" %s %f,%f\n","angular cuts on pions and protons = ", CUTS.pi1cut,CUTS.pi2cut);
122 else
123 printf(" %s %f,%f\n","angular cuts on pions = ", CUTS.pi1cut,CUTS.pi2cut);
124
125 if (FLAGS.tagged == 0) {
126 if (FLAGS.pion == 0)
127 printf(" %s %f %s\n","min. muons-tagged photon inv.mass^2 = ", CUTS.q2min," GeV^2");
128 else if (FLAGS.pion == 4)
129 printf(" %s %f %s\n","min. protons-tagged photon inv.mass^2 = ", CUTS.q2min," GeV^2" );
130 else if (FLAGS.pion == 5)
131 printf(" %s %f %s\n","min. neutrons-tagged photon inv.mass^2 = ", CUTS.q2min," GeV^2" );
132 else if ((FLAGS.pion == 6)||(FLAGS.pion == 7))
133 printf(" %s %f %s\n","min. kaons-tagged photon inv.mass^2 = ", CUTS.q2min," GeV^2" );
134 else if (FLAGS.pion == 9)
135 printf(" %s %f %s\n"," min. lambdas-tagged photon inv.mass^2 = ", CUTS.q2min," GeV^2" );
136 else
137 printf(" %s %f %s\n","min. pions-tagged photon inv.mass^2 = ", CUTS.q2min," GeV^2" );
138 }
139
140// --- book histograms -----------------------------
141
142 INITHISTO();
143
144// --- set cuts ------------------------------------
145 if (FLAGS.tagged == 0) {
146 cos1min = cos(CUTS.phot2cut*CTES.pi/180.0); // photon1 angle cuts in the
147 cos1max = cos(CUTS.phot1cut*CTES.pi/180.0); // LAB rest frame
148 } else {
149 cos1min = -1.0;
150 cos1max = 1.0;
151 CUTS.gmin = CUTS.gmin/2.0;
152 }
153 cos2min = -1.0; // photon2 angle limits
154 cos2max = 1.0; //
155 cos3min = -1.0; // hadrons/muons angle limits
156 cos3max = 1.0; // in their rest frame
157 if (FLAGS.pion == 0) // virtual photon energy cut
158 qqmin = 4.0*CTES.mmu*CTES.mmu;
159 else if (FLAGS.pion == 1)
160 qqmin = 4.0*CTES.mpi*CTES.mpi;
161 else if (FLAGS.pion == 2)
162 qqmin = 4.0*(CTES.mpi+CTES.mpi0)*(CTES.mpi+CTES.mpi0);
163 else if (FLAGS.pion == 3)
164 qqmin = 16.0*CTES.mpi*CTES.mpi;
165 else if (FLAGS.pion == 4)
166 qqmin = 4.0*CTES.mp*CTES.mp;
167 else if (FLAGS.pion == 5)
168 qqmin = 4.0*CTES.mnt*CTES.mnt;
169 else if (FLAGS.pion == 6)
170 qqmin = 4.0*CTES.mKp*CTES.mKp;
171 else if (FLAGS.pion == 7)
172 qqmin = 4.0*CTES.mKn*CTES.mKn;
173 else if (FLAGS.pion == 8)
174 qqmin = (2.0*CTES.mpi+CTES.mpi0)*(2.0*CTES.mpi+CTES.mpi0);
175 else if (FLAGS.pion == 9)
176 qqmin = 4.0*CTES.mlamb*CTES.mlamb;
177 // else
178 // continue;
179
180 qqmax = CTES.Sp-2.0*sqrt(CTES.Sp)*CUTS.gmin; // if only one photon
181 if (CUTS.q2_max_c < qqmax)
182 qqmax=CUTS.q2_max_c; // external cuts
183
184 // -------------------
185 if ( (CUTS.q2_min_c > qqmin) && (CUTS.q2_min_c < (CTES.Sp*(1.0-2.0*(CUTS.gmin/sqrt(CTES.Sp)+CUTS.w))) ) )
186 qqmin = CUTS.q2_min_c;
187 else {
188 cerr << "------------------------------" << endl;
189 cerr << " Q^2_min TOO SMALL" << endl;
190 cerr << " Q^2_min CHANGE BY PHOKHARA = " << qqmin << " GeV^2" << endl;
191 cerr << "------------------------------" << endl;
192 }
193 // -------------------
194 if(qqmax <= qqmin){
195 cerr << " Q^2_max to small " << endl;
196 cerr << " Q^2_max = " << qqmax << endl;
197 cerr << " Q^2_min = " << qqmin << endl;
198 return 0;
199 }
200
201 // -------------------
202 if (FLAGS.pion == 0) {
203 printf(" %s %f %s\n", "minimal muon-pair invariant mass^2 = ", qqmin," GeV^2");
204 printf(" %s %f %s\n", "maximal muon-pair invariant mass^2 = ", qqmax," GeV^2");
205 } else if (FLAGS.pion == 1) {
206 printf(" %s %f %s\n", "minimal pion-pair invariant mass^2 = ", qqmin," GeV^2");
207 printf(" %s %f %s\n", "maximal pion-pair invariant mass^2 = ", qqmax," GeV^2");
208 } else if (FLAGS.pion == 4) {
209 printf(" %s %f %s\n", "minimal proton-pair invariant mass^2 = ", qqmin," GeV^2");
210 printf(" %s %f %s\n", "maximal proton-pair invariant mass^2 = ", qqmax," GeV^2");
211 } else if (FLAGS.pion == 5) {
212 printf(" %s %f %s\n", "minimal neutron-pair invariant mass^2 = ", qqmin," GeV^2");
213 printf(" %s %f %s\n", "maximal neutron-pair invariant mass^2 = ", qqmax," GeV^2");
214 } else if ((FLAGS.pion == 6) || (FLAGS.pion == 7)) {
215 printf(" %s %f %s\n", "minimal kaon-pair invariant mass^2 = ", qqmin," GeV^2");
216 printf(" %s %f %s\n", "maximal kaon-pair invariant mass^2 = ", qqmax," GeV^2");
217 } else if(FLAGS.pion == 8){
218 printf(" %s %f %s\n", "minimal three-pion invariant mass^2 = ", qqmin," GeV^2");
219 printf(" %s %f %s\n", "maximal three-pion invariant mass^2 = ", qqmax," GeV^2");
220 } else if(FLAGS.pion == 9){
221 printf(" %s %f %s\n", "minimal lambda-pair invariant mass^2 = ", qqmin," GeV^2");
222 printf(" %s %f %s\n", "maximal lambda-pair invariant mass^2 = ", qqmax," GeV^2");
223 } else {
224 printf(" %s %f %s\n", "minimal four-pion invariant mass^2 = ", qqmin," GeV^2" );
225 printf(" %s %f %s\n", "maximal four-pion invariant mass^2 = ", qqmax," GeV^2");
226 }
227 // -------------------
228 if (FLAGS.nlo == 0) {
229 cout << "Born" << endl;
230 if(FLAGS.fsrnlo != 0){
231 cerr << "WRONG FSRNLO flag - only fsrnlo=0 allowed for Born" << endl;
232 return 0;
233 }
234 }
235 // -------------------
236 if((FLAGS.pion == 9) && (FLAGS.nlo != 0)) {
237 cerr << "WRONG NLO flag - only Born allowed for Lambdas"<< endl;
238 cerr << "If you feel that you need NLO"<< endl;
239 cerr << "please contact the authors"<< endl;
240 return 0;
241 }
242 // -------------------
243 if (FLAGS.nlo == 1) printf(" %s %f\n", "NLO: soft photon cutoff w = ",CUTS.w);
244 if ((FLAGS.pion <= 1) || (FLAGS.pion == 6))
245 {
246
247 if( ! ((FLAGS.fsr == 1) || (FLAGS.fsr == 2) || (FLAGS.fsrnlo == 0)
248 || (FLAGS.fsr == 1) || (FLAGS.fsrnlo == 1)
249 || (FLAGS.fsr == 0) || (FLAGS.fsrnlo == 0))) {
250 cerr << "WRONG combination of FSR, FSRNLO flags" <<endl;
251 return 0;
252 }
253
254 // ------------------
255 if (FLAGS.fsr == 0)
256 cout << "ISR only" << endl;
257 else if (FLAGS.fsr == 1)
258 cout << "ISR+FSR" << endl;
259 else if (FLAGS.fsr == 2) {
260 if (FLAGS.nlo == 0)
261 cout << "ISR+INT+FSR" << endl;
262 else {
263 cerr << "WRONG FSR flag: interference is included only for nlo=0" << endl;
264 return 0;
265 }
266 }
267 else {
268 cerr << "WRONG FSR flag" << FLAGS.fsr << endl;
269 return 0;
270 }
271
272 if(FLAGS.fsrnlo == 1)
273 cout << "IFSNLO included" << endl;
274 }
275 else
276 {
277 if((FLAGS.fsr == 0) && (FLAGS.fsrnlo == 0))
278 cout << "ISR only" << endl;
279 else {
280 cerr << "FSR is implemented only for pi+pi-, mu+mu- and K+K- modes" << endl;
281 return 0;
282 }
283 }
284
285 // ------------------
286 if(FLAGS.ivac == 0){
287 cout << "Vacuum polarization is NOT included" << endl;
288 } else if(FLAGS.ivac == 1){
289 cout << "Vacuum polarization is included" << endl;
290 } else {
291 cout << "WRONG vacuum polarization switch" << endl;
292 return 0;
293 }
294// -----------------
295 if(FLAGS.pion == 1){
296 if(FLAGS.FF_pion == 0)
297 cout << "Kuhn-Santamaria PionFormFactor" << endl;
298 else if(FLAGS.FF_pion == 1)
299 cout << "Gounaris-Sakurai PionFormFactor" << endl;
300 else {
301 cout << "WRONG PionFormFactor switch" << endl;
302 return 0;
303 }
304// ------
305 if(FLAGS.fsr != 0){
306 if(FLAGS.f0_model == 0)
307 cout << "f0+f0(600): K+K- model" << endl;
308 else if(FLAGS.f0_model == 1)
309 cout << "f0+f0(600): \"no structure\" model" << endl;
310 else if(FLAGS.f0_model == 2)
311 cout << "NO f0+f0(600)" << endl;
312 else if(FLAGS.f0_model == 3)
313 cout << "only f0, KLOE: Cesare Bini-private communication" << endl;
314 else {
315 cout << "WRONG f0+f0(600) switch" << endl;
316 return 0;
317 }
318 }
319 }
320//
321// =================================================
322// --- finding the maximum -------------------------
323 k = nm;
324 for( i = 0; i<2; i++)
325 {
326 MAXIMA.Mmax[i] = 1.0;
327 MAXIMA.gross[i] = 0.0;
328 MAXIMA.klein[i] = 0.0;
329 }
330 if (FLAGS.nlo == 0)
331 MAXIMA.Mmax[1]=0.0; // only 1 photon events generated
332
333 for( i = 1; i<=2; i++) { // initializing the MC loop
334 MAXIMA.tr[0] = 0.0;
335 MAXIMA.tr[1] = 0.0;
336 MAXIMA.count[0] = 0.0;
337 MAXIMA.count[1] = 0.0;
338
339 // =================================================
340 // --- beginning the MC loop event generation ------
341 for(j = 1; j <= k; j++)
342 {
343 RANLXDF(Ar_r,1);
344 Ar[1] = Ar_r[0];
345
346 if (Ar[1] <= (MAXIMA.Mmax[0]/(MAXIMA.Mmax[0]+MAXIMA.Mmax[1]))) {
347 MAXIMA.count[0] = MAXIMA.count[0]+1.0;
348 GEN_1PH(i,qqmin,qqmax,cos1min,cos1max,cos3min,cos3max);
349 }
350 else {
351 MAXIMA.count[1] = MAXIMA.count[1]+1.0;
352 GEN_2PH(i,qqmin,cos1min,cos1max,cos2min,cos2max,cos3min,cos3max);
353 }
354 }
355 // --- end of the MC loop --------------------------
356 // =================================================
357 // --- for the second run ---
358 k = nges;
359 if (i == 1) {
360 MAXIMA.Mmax[0] = MAXIMA.gross[0]+.05*sqrt(MAXIMA.gross[0]*MAXIMA.gross[0]);
361 MAXIMA.Mmax[1] = MAXIMA.gross[1]+(.03+.02*CTES.Sp)*sqrt(MAXIMA.gross[1]*MAXIMA.gross[1]);
362
363 if((FLAGS.pion == 1) && (FLAGS.fsrnlo == 1))
364 MAXIMA.Mmax[1]=MAXIMA.Mmax[1]*1.5;
365 if((FLAGS.pion == 0) && (FLAGS.fsrnlo == 1))
366 MAXIMA.Mmax[1]=MAXIMA.Mmax[1]*1.5;
367
368 if((FLAGS.pion == 0) && (FLAGS.fsr == 1) && (FLAGS.fsrnlo == 0))
369 MAXIMA.Mmax[1]=MAXIMA.Mmax[1]*1.2;
370
371 if((FLAGS.pion == 2) || (FLAGS.pion == 3)){
372 MAXIMA.Mmax[0]=MAXIMA.Mmax[0]*1.1;
373 MAXIMA.Mmax[1]=MAXIMA.Mmax[1]*1.1;
374 }
375
376 if(FLAGS.pion == 8){
377 MAXIMA.Mmax[0]=MAXIMA.Mmax[0]*1.08;
378 MAXIMA.Mmax[1]=MAXIMA.Mmax[1]*1.1;
379 }
380 }
381 }
382 // --- end of the second run -----------------------
383
384 // =================================================
385 if(FLAGS.pion == 9)
386 MAXIMA.Mmax[1] = MAXIMA.Mmax[1] * (1.0 + LAMBDA_PAR.alpha_lamb)*(1.0 + LAMBDA_PAR.alpha_lamb)
387 * LAMBDA_PAR.ratio_lamb*LAMBDA_PAR.ratio_lamb;
388
389 // --- save histograms -----------------------------
390
391 ENDHISTO();
392
393 // --- value of the cross section ------------------
394 if (FLAGS.nlo == 0) {
395 sigma = MAXIMA.Mmax[0]/MAXIMA.count[0]*MAXIMA.tr[0];
396 dsigm = MAXIMA.Mmax[0]*sqrt((MAXIMA.tr[0]/MAXIMA.count[0]-(MAXIMA.tr[0]/MAXIMA.count[0])*(MAXIMA.tr[0]/MAXIMA.count[0]))/MAXIMA.count[0]);
397 } else {
398 sigma1 = MAXIMA.Mmax[0]/MAXIMA.count[0]*MAXIMA.tr[0];
399 sigma2 = MAXIMA.Mmax[1]/MAXIMA.count[1]*MAXIMA.tr[1];
400 dsigm1 = MAXIMA.Mmax[0]*sqrt((MAXIMA.tr[0]/MAXIMA.count[0]-(MAXIMA.tr[0]/MAXIMA.count[0])*(MAXIMA.tr[0]/MAXIMA.count[0]))/MAXIMA.count[0]);
401 dsigm2 = MAXIMA.Mmax[1]*sqrt((MAXIMA.tr[1]/MAXIMA.count[1]-(MAXIMA.tr[1]/MAXIMA.count[1])*(MAXIMA.tr[1]/MAXIMA.count[1]))/MAXIMA.count[1]);
402
403 sigma = sigma1+sigma2;
404 dsigm = sqrt(dsigm1*dsigm1 + dsigm2*dsigm2);
405 }
406
407// --- output --------------------------------------
408 cout << "----------------------------------------------------" << endl;
409 cout << int(MAXIMA.tr[0]+MAXIMA.tr[1]) << " total events accepted of " << endl;
410 cout << int(nges) << " total events generated" << endl;
411 cout << int(MAXIMA.tr[0]) << " one photon events accepted of " << endl;
412 cout << int(MAXIMA.count[0]) << " events generated" << endl;
413 cout << int(MAXIMA.tr[1]) << " two photon events accepted of " << endl;
414 cout << int(MAXIMA.count[1]) << " events generated" << endl;
415 cout << endl;
416 if (FLAGS.nlo != 0) cout << "sigma1(nbarn) = " << sigma1 << " +- " << dsigm1 << endl;
417 if (FLAGS.nlo != 0) cout << "sigma2(nbarn) = " << sigma2 << " +- " << dsigm2 << endl;
418 cout << "sigma (nbarn) = " << sigma << " +- " <<dsigm << endl;
419 cout << endl;
420 cout << "maximum1 = " << MAXIMA.gross[0] << " minimum1 = " << MAXIMA.klein[0] << endl;
421 cout << "Mmax1 = " << MAXIMA.Mmax[0] << endl;
422 cout << "maximum2 = " << MAXIMA.gross[1] << " minimum2 = " << MAXIMA.klein[1] << endl;
423 cout << "Mmax2 = " << MAXIMA.Mmax[1] << endl;
424 cout << "----------------------------------------------------" << endl;
425
426// --- saves the new seed --------------------------
427/*
428 close (9,DISP="delete")
429 open(9,file="seed.dat",type="new")
430
431 call rlxdgetf(s_seed)
432 write(9,*)s_seed
433
434 end
435 */
436
437return 0;
438}
cos
double cos(const BesAngle a)
Definition BesAngle.h:213
FLAGS
#define FLAGS
Definition BesBdkRc.cxx:82
sigma
TTree * sigma
Definition CurConst_dEdXCal.cxx:6
CTES
#define CTES
Definition PhokharaDef.h:20
LAMBDA_PAR
#define LAMBDA_PAR
Definition PhokharaDef.h:61
MAXIMA
#define MAXIMA
Definition PhokharaDef.h:84
CUTS
#define CUTS
Definition PhokharaDef.h:30

◆ PROTOCCALLSFSUB0() [1/3]

PROTOCCALLSFSUB0 ( ENDHISTO ,
endhisto  )

◆ PROTOCCALLSFSUB0() [2/3]

PROTOCCALLSFSUB0 ( INITHISTO ,
inithisto  )

◆ PROTOCCALLSFSUB0() [3/3]

PROTOCCALLSFSUB0 ( WRITEEVENT ,
writeevent  )

◆ PROTOCCALLSFSUB1()

PROTOCCALLSFSUB1 ( RLXDRESETF ,
rlxdresetf ,
INTV  )

◆ PROTOCCALLSFSUB2()

PROTOCCALLSFSUB2 ( RANLXDF ,
ranlxdf ,
DOUBLEV ,
INT  )

◆ PROTOCCALLSFSUB3()

PROTOCCALLSFSUB3 ( INPUT ,
input ,
PLONG ,
PINT ,
PSTRING  )

◆ PROTOCCALLSFSUB7()

PROTOCCALLSFSUB7 ( GEN_1PH ,
gen_1ph_ ,
INT ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE  )

◆ PROTOCCALLSFSUB8()

PROTOCCALLSFSUB8 ( GEN_2PH ,
gen_2ph_ ,
INT ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE ,
PDOUBLE  )