Index: /trunk/WuerzburgSoft/Thomas/mphys/Changelog =================================================================== --- /trunk/WuerzburgSoft/Thomas/mphys/Changelog (revision 1366) +++ /trunk/WuerzburgSoft/Thomas/mphys/Changelog (revision 1367) @@ -8,4 +8,6 @@ - changed the eps of all integrals to 1e-2 - changed the p_e function to use the Compton function + - made all function which are used for integration only static + and changed the function definition from=0, to=0 * MPhoton.[h,cc]: @@ -25,4 +27,6 @@ * phys.C: - changed the integral eps to 1e-2 + - small changed to the output + - removed integral==0, replaced by Ep==pow(10,0) Index: /trunk/WuerzburgSoft/Thomas/mphys/MElectron.cc =================================================================== --- /trunk/WuerzburgSoft/Thomas/mphys/MElectron.cc (revision 1366) +++ /trunk/WuerzburgSoft/Thomas/mphys/MElectron.cc (revision 1367) @@ -90,5 +90,6 @@ return DiSum(x); - TF1 IntLi("Li", Li, 0, z, 0); + // TF1 IntLi("Li", Li, 0, z, 0); + static TF1 IntLi("Li", Li, 0, 0, 0); const Double_t integ = IntLi.Integral(0, z, (Double_t*)NULL, 1e-2); return -integ; @@ -160,9 +161,4 @@ const Double_t z = k ? k[0] : 0; - Double_t val[3] = { E[0], z }; // E[GeV] - - const Double_t from = 1e-17; - const Double_t to = 1e-11; - /* -------------- old ---------------- Double_t from = 1e-15; @@ -171,5 +167,10 @@ ----------------------------------- */ - TF1 func("Compton", Compton, from, to, 2); // [0, inf] + static TF1 func("Compton", Compton, 0, 0, 2); // [0, inf] + + const Double_t from = 1e-17; + const Double_t to = 1e-11; + + Double_t val[3] = { E[0], z }; // E[GeV] Double_t integ = func.Integral(from, to, val, 1e-2); // [Gev] [0, inf] @@ -260,5 +261,4 @@ const Double_t Gamma = 4.*omega; - // --old-- fQ.SetRange(1e-6, 1./(1.+ 2.*Gamma)); fQ.SetParameter(0, Gamma); Index: /trunk/WuerzburgSoft/Thomas/mphys/MPairProduction.cc =================================================================== --- /trunk/WuerzburgSoft/Thomas/mphys/MPairProduction.cc (revision 1366) +++ /trunk/WuerzburgSoft/Thomas/mphys/MPairProduction.cc (revision 1367) @@ -16,8 +16,7 @@ ! ! -! Author(s): Harald Kornmayer 1/2001 (harald@mppmu.mpg.de) -! Author(s): Thomas Bretz 12/2000 (tbretz@uni-sw.gwdg.de) +! Author(s): Thomas Bretz 12/2000 ! -! Copyright: MAGIC Software Development, 2000-2001 +! Copyright: MAGIC Software Development, 2000-2002 ! ! Index: /trunk/WuerzburgSoft/Thomas/mphys/MPhoton.cc =================================================================== --- /trunk/WuerzburgSoft/Thomas/mphys/MPhoton.cc (revision 1366) +++ /trunk/WuerzburgSoft/Thomas/mphys/MPhoton.cc (revision 1367) @@ -92,11 +92,10 @@ Double_t val[2] = { Eg, Ep }; - const Double_t from = -1.0; - const Double_t to = 1.-E0*E0/(2.*Eg*Ep); // Originally Was: 1. + static TF1 f("int1", Int1, 0, 0, 2); - TF1 f("int1", Int1, from, to, 2); - + const Double_t from = -1.0; + const Double_t to = 1.-E0*E0/(2.*Eg*Ep); // Originally Was: 1. const Double_t int1 = f.Integral(from, to, val, 1e-2); // [m^2] - const Double_t planck = MParticle::Planck(&Ep, &z); // [GeV^2] + const Double_t planck = MParticle::Planck(&Ep, &z); // [GeV^2] const Double_t res = planck * int1; @@ -124,10 +123,9 @@ Double_t val[2] = { Eg, z }; + static TF1 f("int2", Int2, 0, 0, 2); + Double_t lolim = E0*E0/Eg; Double_t inf = Eg<1e6 ? 3e-11*(z+1) : 3e-12*(z+1); - - TF1 f("int2", Int2, lolim, inf, 2); - - Double_t int2 = f.Integral(lolim, inf, val, 1e-2); //[GeV^3 m^2] + Double_t int2 = f.Integral(lolim, inf, val, 1e-2); //[GeV^3 m^2] if (int2==0) Index: /trunk/WuerzburgSoft/Thomas/mphys/phys.C =================================================================== --- /trunk/WuerzburgSoft/Thomas/mphys/phys.C (revision 1366) +++ /trunk/WuerzburgSoft/Thomas/mphys/phys.C (revision 1367) @@ -61,4 +61,6 @@ TF1 func("RndTheta", Sbar_sigmas, 0, log10(f), 0); + func.SetNpx(50); + Double_t sbar = pow(10, func.GetRandom()); Double_t theta = acos(1.-sbar*2/f); @@ -116,6 +118,6 @@ void DoIt() { - Double_t startz = 0.003; //MParticle::ZofR(&R); - Double_t R = MParticle::RofZ(&startz); // [kpc] + Double_t R = 100; //MParticle::RofZ(&startz); // [kpc] + Double_t startz = MParticle::ZofR(&R); cout << "R = " << R << "kpc" << endl; @@ -127,5 +129,5 @@ Double_t runtime = 15*60; // [s] - Double_t lo = 1e4; + Double_t lo = 1e4; Double_t hi = 1e11; Double_t alpha = -2; @@ -274,4 +276,5 @@ phot.SetParameter(0, Eg); phot.SetParameter(1, p->GetZ()); + /* if (phot.Integral(-log10(Eg)-8, -10.5, (Double_t*)NULL, 1e-2)==0) { @@ -279,6 +282,11 @@ continue; } - + */ Double_t Ep = pow(10, phot.GetRandom()); + if (Ep==pow(10,0)) + { + cout << "z" << flush; + continue; + } Double_t theta = RandomTheta(Eg, Ep); if (theta==0) @@ -339,5 +347,8 @@ else { - cout << "i" << flush; // ignored + if (p->GetEnergy()<=lo) + cout << "e" << flush; + else + cout << "t" << flush; // ignored delete p; } @@ -374,5 +385,5 @@ cout << "Created " << n << " gammas from source with E^" << alpha << endl; - cout << "Processed " << Form("%.1f", n/(TStopwatch::GetRealTime()-starttime)) << " gammas/sec." << endl; + cout << "Processing time: " << Form("%.1f", (TStopwatch::GetRealTime()-starttime)/n) << " sec/gamma." << endl; // ------------------------------