#include <Lfluk.h>
Public Member Functions | |
| Lfluk () | |
| virtual | ~Lfluk () |
| void | SetConfig (const std::string param, const std::string skew, const std::string na49fit, const double wmin, const double wmax) |
| void | Config (const Registry ®) |
| void | Init () |
| void | SetParameters (std::vector< double > par) |
| const std::vector< double > & | GetParameters () |
| double | GetWeight (int ptype, double pt, double pz) |
| double | GetMeanPT (int ptype) |
| double | GetReweightedMeanPT (int ptype) |
| double | GetPTshift (int ptype, double pz_low=0., double pz_up=120.) |
| double | GetNFrac (int ptype, double pz_low=0., double pz_up=120., double pt_low=0., double pt_up=1.) |
| TH2F * | GetReweightedPTXF (int ptype) |
| TH2F * | GetPTXF (int ptype) |
| TH2D * | GetWeightHistogram (int ptype) |
| TH2 * | GetWeightHist (const std::string &ptype) |
| TH1 * | GetFlukKPiRatio () |
| TH1 * | GetPiPlusReWProj () |
| TH1 * | GetKPlusf05Proj () |
| TH1 * | GetKPlusf05PipReWRatio () |
Private Types | |
| typedef enum Lfluk::EParticleType | ParticleType_t |
| enum | EParticleType { kPiPlus = 8, kPiMinus = 9, kKPlus = 11, kKMinus = 12, kK0L = 10, kUnknown = 0 } |
Private Member Functions | |
| double | CalcWeight (int ptype, double pt, double pz) |
| double | GetSKZPWeight (const Lfluk::ParticleType_t particle, const double pT, const double xF) |
| double | GetSKPCWeight (const Lfluk::ParticleType_t particle, const double pT, const double xF) |
| double | GetQuadWeight (const Lfluk::ParticleType_t particle, const double pT, const double xF) |
| double | GetMeanPTWeight (const Lfluk::ParticleType_t particle, const double pT, const double xF) |
| void | SkewSKZP (const double xF, const double A, const double B, const double C, long double &Ap, long double &Bp, long double &Cp) |
| void | SkewSKPC (const double xF, const double A, const double B, const double C, const double D, const double E, double &Ap, double &Bp, double &Cp, double &Dp, double &Ep) |
| double | GetSKZPA (const Lfluk::ParticleType_t particle, double xf) const |
| double | GetSKZPB (const Lfluk::ParticleType_t particle, double xf) const |
| double | GetSKZPC (const Lfluk::ParticleType_t particle, double xf) const |
| double | GetSKPCA (const Lfluk::ParticleType_t particle, double xf) const |
| double | GetSKPCB (const Lfluk::ParticleType_t particle, double xf) const |
| double | GetSKPCC (const Lfluk::ParticleType_t particle, double xf) const |
| double | GetSKPCD (const Lfluk::ParticleType_t particle, double xf) const |
| double | GetSKPCE (const Lfluk::ParticleType_t particle, double xf) const |
| std::string | GetParticleName (Lfluk::ParticleType_t ptype) |
| Lfluk::ParticleType_t | GetParticleEnum (const std::string &ptype) |
| Lfluk::ParticleType_t | GetParticleEnum (int ptype) |
| void | RecalculateWeights () |
Private Attributes | |
| std::string | fParam |
| std::string | fSkew |
| std::string | fNA49Fit |
| std::string | fPath |
| double | fMinW |
| double | fMaxW |
| double | fKPiRatScale |
| bool | fWghtFromHist |
| bool | fUseFlukPi |
| std::vector< double > | fPar |
| std::vector< ParticleType_t > | fPlist |
| std::map< Lfluk::ParticleType_t, TH2F * > | fPTXF |
| std::map< Lfluk::ParticleType_t, TH2F * > | fWeightedPTXF |
| std::map< Lfluk::ParticleType_t, TH2D * > | fWeightHist |
| std::map< Lfluk::ParticleType_t, double > | fMeanPT |
| std::map< Lfluk::ParticleType_t, double > | fN |
| std::map< Lfluk::ParticleType_t, double > | fNWeighted |
| std::map< Lfluk::ParticleType_t, double > | fWeightedMeanPT |
| std::map< Lfluk::ParticleType_t, int > | fNBinsY |
| std::map< Lfluk::ParticleType_t, int > | fNBinsX |
| TF1 * | fNA49pirat |
| std::vector< double > | fPiRatio |
| std::vector< double > | fFluka05KPiRatio |
| std::vector< double > | fFluka05PiRatio |
| std::vector< double > | fKPiRatio |
| bool | isFirst |
| TH1 * | fFlukKPiRat |
| TH1 * | fFlukPiRat |
| TH1 * | fPiPlusReWProj |
| TH1 * | fKPlusf05Proj |
| TH1 * | fKPlusf05PipReWRatio |
|
|
Referenced by CalcWeight(), GetNFrac(), GetParticleEnum(), GetPTshift(), GetWeight(), and RecalculateWeights(). |
|
|
Definition at line 56 of file Lfluk.h. 00057 {
00058 kPiPlus = 8,
00059 kPiMinus = 9,
00060 kKPlus = 11,
00061 kKMinus = 12,
00062 kK0L = 10,
00063 kUnknown = 0
00064 } ParticleType_t;
|
|
|
Definition at line 21 of file Lfluk.cxx. 00022 :fParam("skzp"), 00023 fSkew("lin"), 00024 fNA49Fit("Poly3"), 00025 fPath(""), 00026 fMinW(-999.0), 00027 fMaxW(999.0), 00028 fKPiRatScale(1.0), 00029 fWghtFromHist(true), 00030 fUseFlukPi(false), 00031 fPar(0), 00032 fFlukKPiRat(0) 00033 { 00034 00035 }
|
|
|
Definition at line 38 of file Lfluk.cxx. 00039 {
00040 }
|
|
||||||||||||||||
|
Definition at line 528 of file Lfluk.cxx. References fKPiRatio, fN, fNWeighted, fPar, fParam, fPiRatio, GetMeanPTWeight(), GetParticleEnum(), GetQuadWeight(), GetSKPCWeight(), GetSKZPWeight(), kKPlus, kPiPlus, MAXMSG, and ParticleType_t. Referenced by GetWeight(), and RecalculateWeights(). 00529 {
00530 double weight=1.;
00531
00532 if (fPar.size()==0)
00533 {
00534 MAXMSG("MCReweight",Msg::kWarning,10)
00535 <<"You need to set the parameters before calling "
00536 <<"Zfluk::CalcWeight (use SetParameters(vector<double>))"<<std::endl
00537 <<"Returning weight = "<<weight<<std::endl;
00538 return weight;
00539 }
00540
00541 double xF=pz/120.;
00542 Lfluk::ParticleType_t eptype=GetParticleEnum(ptype);
00543
00544 // This is SKZP parameterization
00545 if (xF>1.||xF<0.) return weight;
00546
00547 if (eptype==kPiPlus)
00548 {
00549 if (pT>1.||pT<0.) return weight;
00550
00551 if(fParam.find("skzp") != string::npos) //default
00552 {
00553 weight = Lfluk::GetSKZPWeight(eptype, pT, xF);
00554 }
00555 if(fParam.find("skpc") != string::npos)
00556 {
00557 weight = Lfluk::GetSKPCWeight(eptype, pT, xF);
00558 }
00559 if(fParam.find("quad") != string::npos)
00560 {
00561 weight = Lfluk::GetQuadWeight(eptype, pT, xF);
00562 }
00563 if(fParam.find("meanpt") != string::npos)
00564 {
00565 weight = Lfluk::GetMeanPTWeight(eptype, pT, xF);
00566 }
00567
00568 if(weight < fMinW) weight = fMinW;
00569 if(weight > fMaxW) weight = fMaxW;
00570
00571 }
00572 else if (eptype==kKPlus)
00573 {
00574 weight=fKPiRatio[int(pz)+1];
00575 }
00576 else if (eptype==kPiMinus)
00577 {
00578 weight=CalcWeight(kPiPlus,pT,pz)*fPiRatio[int(pz)+1];
00579 }
00580 else if (eptype==kKMinus)
00581 {
00582 weight=CalcWeight(kKPlus,pT,pz);
00583 }
00584 else if (eptype==kK0L)
00585 {
00586 //N(K0L) approximately given with (N(K+)+3*N(K-))/4
00587 weight= ((fNWeighted[kKPlus]+3.*fNWeighted[kKMinus])
00588 /(fN[kKPlus]+3.*fN[kKMinus]));
00589 }
00590
00591 if(weight < -999.0) weight = fMinW;
00592 if(weight > 999.0) weight = fMaxW;
00593
00594
00595 return weight;
00596 }
|
|
|
Definition at line 273 of file Lfluk.cxx. References fKPiRatScale, fMaxW, fMinW, fNA49Fit, fParam, fPath, fSkew, fUseFlukPi, fWghtFromHist, Registry::Get(), and reg. 00274 {
00275 reg.Get("WEIGHT_MIN", fMinW);
00276 reg.Get("WEIGHT_MAX", fMaxW);
00277 reg.Get("KPiRatioScale", fKPiRatScale);
00278
00279 if(fMaxW < 990.0 && fMinW > -990.0 && fMaxW <= fMinW)
00280 {
00281 fMaxW = fMinW + 10.0;
00282 cout << "Lfluk::Config - Weight Max <= Weight Min, setting fMaxW = " << fMaxW << " fMinW = " << fMinW << endl;
00283 }
00284
00285 const char *value_char = 0;
00286
00287 if(reg.Get("ZRatTuneParam", value_char) && value_char)
00288 {
00289 stringstream strstrm;
00290 strstrm << value_char;
00291 fParam = strstrm.str();
00292 }
00293
00294 value_char = 0;
00295 if(reg.Get("ZRatTuneSkew", value_char) && value_char)
00296 {
00297 stringstream strstrm;
00298 strstrm << value_char;
00299 fSkew = strstrm.str();
00300 }
00301
00302 value_char = 0;
00303 if(reg.Get("NA49FIT", value_char) && value_char)
00304 {
00305 stringstream strstrm;
00306 strstrm << value_char;
00307 fNA49Fit = strstrm.str();
00308 }
00309
00310 value_char = 0;
00311 if(reg.Get("ROOTFILES", value_char) && value_char)
00312 {
00313 stringstream strstrm;
00314 strstrm << value_char;
00315 fPath = strstrm.str();
00316 }
00317
00318 value_char = 0;
00319 if(reg.Get("WGHTFROMHIST", value_char) && value_char)
00320 {
00321 if(std::strcmp(value_char, "no") == 0)
00322 {
00323 fWghtFromHist = false;
00324 }
00325 }
00326
00327 value_char = 0;
00328 if(reg.Get("UseFlukaPiRatio", value_char) && value_char)
00329 {
00330 if(std::strcmp(value_char, "yes") == 0)
00331 {
00332 fUseFlukPi = true;
00333 }
00334 }
00335
00336
00337 cout << "Lfluk::Config" << endl;
00338 cout << " Root File Path = " << fPath << endl
00339 << " Get Weight From Hist = " << fWghtFromHist << endl
00340 << " Parameterization = " << fParam << endl
00341 << " Skew Config = " << fSkew << endl
00342 << " Use Fluka Pi Ratio = " << fUseFlukPi << endl
00343 << " NA49 Ratio Fit = " << fNA49Fit << endl
00344 << " KPi Ratio Scale = " << fKPiRatScale << endl
00345 << " Weight Min = " << fMinW << endl
00346 << " Weight Max = " << fMaxW << endl;
00347
00348
00349 }
|
|
|
Definition at line 49 of file Lfluk.h. 00049 {return fFlukKPiRat;};
|
|
|
Definition at line 52 of file Lfluk.h. 00052 {return fKPlusf05PipReWRatio;};
|
|
|
Definition at line 51 of file Lfluk.h. 00051 {return fKPlusf05Proj;};
|
|
|
Definition at line 39 of file Lfluk.h. 00039 {return fMeanPT[GetParticleEnum(ptype)];};
|
|
||||||||||||||||
|
Definition at line 905 of file Lfluk.cxx. References fPar. Referenced by CalcWeight(). 00906 {
00907 double weight = 1.0;
00908
00909 if(particle == kPiPlus)
00910 {
00911 if(fPar.size() != 5)
00912 {
00913 cout << "Problem in Lfluk::GetMeanPTWeight - Number of pars is not 4: npars = " << fPar.size() << endl;
00914 }
00915 else
00916 {
00917 double linxF = fPar[3]+fPar[4]*xF;
00918 if(linxF == 0.0)
00919 {
00920 linxF = 1e-20;
00921 }
00922 double exparg = fPar[2]*pT / linxF;
00923 if(exparg > 100.0) exparg = 100.0;
00924 if(exparg < -100.0) exparg = -100.0;
00925 weight = (fPar[0] + fPar[1]*xF) * exp(exparg);
00926 }
00927 }
00928
00929
00930 //cout << "xF = " << xF << " pT = " << pT << " meanPT weight = " << weight << endl;
00931 return weight;
00932 }
|
|
||||||||||||||||||||||||
|
Definition at line 1330 of file Lfluk.cxx. References fPTXF, fWeightedPTXF, GetParticleEnum(), and ParticleType_t. 01331 {
01332 //returns the fractional change in number of particles due to reweighting in xf-pt region
01333 //given by pz_low, pz_up and pt_low, pt_up
01334
01335 Lfluk::ParticleType_t eptype=GetParticleEnum(ptype);
01336 if (fPTXF.find(eptype)==fPTXF.end()) return 1.;
01337
01338 int xlow=int(pz_low+1.);
01339 int xup=int(pz_up+1.);
01340 int ylow=int(50.*pt_low+1.);
01341 int yup=int(50.*pt_up+1.);
01342
01343 double wsum(0.);
01344 double sum(0.);
01345 for (int i=xlow;i<xup;i++) {
01346 for (int j=ylow;j<yup;j++) {
01347 sum+=fPTXF[eptype]->GetBinContent(i,j);
01348 wsum+=fWeightedPTXF[eptype]->GetBinContent(i,j);
01349 }
01350 }
01351
01352 return wsum/sum;
01353 }
|
|
|
Definition at line 35 of file Lfluk.h. 00035 {return fPar;};
|
|
|
Definition at line 1383 of file Lfluk.cxx. References ParticleType_t. 01384 {
01385 switch (ptype)
01386 {
01387 case 8 : return kPiPlus; break;
01388 case 211 : return kPiPlus; break;
01389 case 9 : return kPiMinus; break;
01390 case -211: return kPiMinus; break;
01391 case 11 : return kKPlus; break;
01392 case 321 : return kKPlus; break;
01393 case 12 : return kKMinus; break;
01394 case -321: return kKMinus; break;
01395 case 10 : return kK0L; break;
01396 case 130 : return kK0L; break;
01397 default: return kUnknown; break;
01398 }
01399 return kUnknown;
01400
01401 }
|
|
|
Definition at line 1373 of file Lfluk.cxx. References ParticleType_t. Referenced by CalcWeight(), GetNFrac(), GetPTshift(), GetWeight(), and GetWeightHist(). 01374 {
01375 if (ptype == "PiPlus") return kPiPlus;
01376 else if(ptype == "PiMinus") return kPiMinus;
01377 else if(ptype == "KPlus") return kKPlus;
01378 else if(ptype == "KMinus") return kKMinus;
01379 else return kUnknown;
01380 }
|
|
|
Definition at line 1356 of file Lfluk.cxx. References kK0L, kKMinus, kKPlus, kPiMinus, kPiPlus, and kUnknown. Referenced by GetWeight(), Init(), and RecalculateWeights(). 01357 {
01358 switch (ptype)
01359 {
01360 case kPiPlus: return "PiPlus"; break;
01361 case kPiMinus: return "PiMinus"; break;
01362 case kKPlus: return "KPlus"; break;
01363 case kKMinus: return "KMinus"; break;
01364 case kK0L: return "K0L"; break;
01365 case kUnknown: return "Unknown"; break;
01366 default: return "Unknown"; break;
01367 }
01368
01369 return "Unknown";
01370 }
|
|
|
Definition at line 50 of file Lfluk.h. 00050 {return fPiPlusReWProj;};
|
|
||||||||||||||||
|
Definition at line 1299 of file Lfluk.cxx. References fMeanPT, fPTXF, fWeightedMeanPT, fWeightedPTXF, GetParticleEnum(), and ParticleType_t. 01300 {
01301 //returns the pT shift for ptype particle due to reweighting in xf slice
01302 //given by pz_low and pz_up
01303 //this function returns correct ptshift only if the parameters were set using SetParameters
01304
01305 Lfluk::ParticleType_t eptype=GetParticleEnum(ptype);
01306 if (fPTXF.find(eptype)==fPTXF.end()) return 0.;
01307
01308 if (pz_low==0.&&pz_up==120.) return fWeightedMeanPT[eptype]-fMeanPT[eptype];
01309
01310 int xlow=int(pz_low+1.);
01311 int xup=int(pz_up+1.);
01312
01313 double pt;
01314 double meanpt(0.);
01315 double sum(0.);
01316 for (int i=xlow;i<xup;i++) {
01317 for (int j=1;j<51;j++) {
01318 pt=0.02*j-0.01;
01319 meanpt+=fWeightedPTXF[eptype]->GetBinContent(i,j)*pt;
01320 sum+=fWeightedPTXF[eptype]->GetBinContent(i,j);
01321 }
01322 }
01323 meanpt/=sum;
01324 meanpt*=1000.; //GeV to MeV
01325
01326 return meanpt-fMeanPT[eptype];
01327
01328 }
|
|
|
Definition at line 44 of file Lfluk.h. Referenced by Init(). 00044 {return fPTXF[GetParticleEnum(ptype)];};
|
|
||||||||||||||||
|
Definition at line 692 of file Lfluk.cxx. Referenced by CalcWeight(). 00693 {
00694 double weight = 1.0;
00695
00696 if(particle == kPiPlus)
00697 {
00698 if(fParam == "quad0")
00699 {
00700 if(fPar.size() != 9)
00701 {
00702 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 9: npars = " << fPar.size() << endl;
00703 }
00704 else
00705 {
00706 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00707 + fPar[3]*pT + fPar[4]*pT*xF + fPar[5]*pT*xF*xF
00708 + fPar[6]*pT*pT + fPar[7]*pT*pT*xF + fPar[8]*pT*pT*xF*xF;
00709 }
00710 }
00711 else if(fParam == "quad1")
00712 {
00713 if(fPar.size() != 8)
00714 {
00715 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00716 }
00717 else
00718 {
00719 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00720 + fPar[3]*pT + fPar[4]*pT*xF + fPar[5]*pT*xF*xF
00721 + fPar[6]*pT*pT + fPar[7]*pT*pT*xF ;
00722 }
00723 }
00724 else if(fParam == "quad2")
00725 {
00726 if(fPar.size() != 8)
00727 {
00728 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00729 }
00730 else
00731 {
00732 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00733 + fPar[3]*pT + fPar[4]*pT*xF + fPar[5]*pT*xF*xF
00734 + fPar[6]*pT*pT + fPar[7]*pT*pT*xF*xF;
00735 }
00736 }
00737 else if(fParam == "quad3")
00738 {
00739 if(fPar.size() != 8)
00740 {
00741 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00742 }
00743 else
00744 {
00745 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00746 + fPar[3]*pT + fPar[4]*pT*xF + fPar[5]*pT*xF*xF
00747 + fPar[6]*pT*pT*xF + fPar[7]*pT*pT*xF*xF;
00748 }
00749 }
00750 else if(fParam == "quad4")
00751 {
00752 if(fPar.size() != 8)
00753 {
00754 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00755 }
00756 else
00757 {
00758 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00759 + fPar[3]*pT + fPar[4]*pT*xF
00760 + fPar[5]*pT*pT + fPar[6]*pT*pT*xF + fPar[7]*pT*pT*xF*xF;
00761 }
00762 }
00763 else if(fParam == "quad5")
00764 {
00765 if(fPar.size() != 8)
00766 {
00767 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00768 }
00769 else
00770 {
00771 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00772 + fPar[3]*pT + fPar[4]*pT*xF*xF
00773 + fPar[5]*pT*pT + fPar[6]*pT*pT*xF + fPar[7]*pT*pT*xF*xF;
00774 }
00775 }
00776 else if(fParam == "quad6")
00777 {
00778 if(fPar.size() != 8)
00779 {
00780 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00781 }
00782 else
00783 {
00784 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00785 + fPar[3]*pT*xF + fPar[4]*pT*xF*xF
00786 + fPar[5]*pT*pT + fPar[6]*pT*pT*xF + fPar[7]*pT*pT*xF*xF;
00787 }
00788 }
00789 else if(fParam == "quad7")
00790 {
00791 if(fPar.size() != 8)
00792 {
00793 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00794 }
00795 else
00796 {
00797 weight = fPar[0] + fPar[1]*xF
00798 + fPar[2]*pT + fPar[3]*pT*xF + fPar[4]*pT*xF*xF
00799 + fPar[5]*pT*pT + fPar[6]*pT*pT*xF + fPar[7]*pT*pT*xF*xF;
00800 }
00801 }
00802 else if(fParam == "quad8")
00803 {
00804 if(fPar.size() != 8)
00805 {
00806 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00807 }
00808 else
00809 {
00810 weight = fPar[0] + fPar[1]*xF*xF
00811 + fPar[2]*pT + fPar[3]*pT*xF + fPar[4]*pT*xF*xF
00812 + fPar[5]*pT*pT + fPar[6]*pT*pT*xF + fPar[7]*pT*pT*xF*xF;
00813 }
00814 }
00815 else if(fParam == "quad9")
00816 {
00817 if(fPar.size() != 8)
00818 {
00819 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 8: npars = " << fPar.size() << endl;
00820 }
00821 else
00822 {
00823 weight = fPar[0]*xF + fPar[1]*xF*xF
00824 + fPar[2]*pT + fPar[3]*pT*xF + fPar[4]*pT*xF*xF
00825 + fPar[5]*pT*pT + fPar[6]*pT*pT*xF + fPar[7]*pT*pT*xF*xF;
00826 }
00827 }
00828 else if(fParam == "quad124")
00829 {
00830 if(fPar.size() != 6)
00831 {
00832 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 7: npars = " << fPar.size() << endl;
00833 }
00834 else
00835 {
00836 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00837 + fPar[3]*pT + fPar[4]*pT*xF
00838 + fPar[5]*pT*pT ;
00839 }
00840 }
00841 else if(fParam == "quad134")
00842 {
00843 if(fPar.size() != 6)
00844 {
00845 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 7: npars = " << fPar.size() << endl;
00846 }
00847 else
00848 {
00849 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00850 + fPar[3]*pT + fPar[4]*pT*xF
00851 + fPar[5]*pT*pT*xF ;
00852 }
00853 }
00854 else if(fParam == "quad1234")
00855 {
00856 if(fPar.size() != 5)
00857 {
00858 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 7: npars = " << fPar.size() << endl;
00859 }
00860 else
00861 {
00862 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00863 + fPar[3]*pT + fPar[4]*pT*xF ;
00864 }
00865 }
00866 else if(fParam == "quad14")
00867 {
00868 if(fPar.size() != 7)
00869 {
00870 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 7: npars = " << fPar.size() << endl;
00871 }
00872 else
00873 {
00874 weight = fPar[0] + fPar[1]*xF + fPar[2]*xF*xF
00875 + fPar[3]*pT + fPar[4]*pT*xF
00876 + fPar[5]*pT*pT + fPar[6]*pT*pT*xF ;
00877 }
00878 }
00879 else if(fParam == "quad147")
00880 {
00881 if(fPar.size() != 6)
00882 {
00883 cerr << "Problem in Lfluk::GetQuadWeight - Number of pars is not 6: npars = " << fPar.size() << endl;
00884 }
00885 else
00886 {
00887 weight = fPar[0] + fPar[1]*xF
00888 + fPar[2]*pT + fPar[3]*pT*xF
00889 + fPar[4]*pT*pT + fPar[5]*pT*pT*xF ;
00890 }
00891 }
00892 else
00893 {
00894 cerr << "Problem in Lfluk::GetQuadWeight - Not a valid Parameterization " << fParam << endl;
00895 }
00896
00897
00898 }
00899
00900 //cout << "xF = " << xF << " pT = " << pT << " quad weight = " << weight << endl;
00901 return weight;
00902 }
|
|
|
Definition at line 40 of file Lfluk.h. 00040 {return fWeightedMeanPT[GetParticleEnum(ptype)];};
|
|
|
Definition at line 43 of file Lfluk.h. Referenced by RecalculateWeights(). 00043 {return fWeightedPTXF[GetParticleEnum(ptype)];};
|
|
||||||||||||
|
Definition at line 1411 of file Lfluk.cxx. References kKMinus, kKPlus, kPiMinus, and kPiPlus. Referenced by GetSKPCWeight(). 01412 {
01413 if (xf<0.25)
01414 {
01415 //
01416 //form: Par[0]*std::pow(1.0-xf, Par[1])*(1.0+Par[2]*xf)*std::pow(xf, -Par[3]);
01417 //
01418 switch (particle)
01419 {
01420 case kPiPlus:
01421 return 0.0462498*std::pow(1.0-xf, 20.8206)*(1.0+11911.3*xf)*std::pow(xf, -1.0*0.96808);
01422 case kPiMinus:
01423 break;
01424 case kKPlus:
01425 break;
01426 case kKMinus:
01427 break;
01428 default:
01429 break;
01430 }
01431 }
01432 else
01433 {
01434 //
01435 //form: Par[4]*(1.0-xf)*std::exp(Par[5]*std::pow(xf, Par[6]));
01436 //
01437 switch (particle)
01438 {
01439 case kPiPlus:
01440 return 1250.0*(1.0-xf)*std::exp(-13.6097*std::pow(xf, 0.597619));
01441 case kPiMinus:
01442 break;
01443 case kKPlus:
01444 break;
01445 case kKMinus:
01446 break;
01447 default:
01448 break;
01449 }
01450 }
01451
01452 return 1.0;
01453
01454 }
|
|
||||||||||||
|
Definition at line 1457 of file Lfluk.cxx. References kKMinus, kKPlus, kPiMinus, and kPiPlus. Referenced by GetSKPCWeight(). 01458 {
01459
01460 //
01461 //form: Par[0]*std::pow(1.0-xf, Par[1])*(1.0+Par[2]*xf)*std::pow(xf, -Par[3]);
01462 //
01463
01464 switch (particle)
01465 {
01466 case kPiPlus:
01467 return 3.31629*std::pow(1.0-xf, 2.48421)*(1.0+22.7285*xf)*std::pow(xf, -1.0*2.38778);
01468 case kPiMinus:
01469 break;
01470 case kKPlus:
01471 break;
01472 case kKMinus:
01473 break;
01474 default:
01475 break;
01476 }
01477
01478 return 0.0;
01479 }
|
|
||||||||||||
|
Definition at line 1482 of file Lfluk.cxx. References kKMinus, kKPlus, kPiMinus, and kPiPlus. Referenced by GetSKPCWeight(). 01483 {
01484 if (xf<0.35)
01485 {
01486 //
01487 //form: Par[0]*std::pow(1.0-xf, Par[1])*(1.0+Par[2]*xf)*std::pow(xf, -Par[3]);
01488 //
01489 switch (particle)
01490 {
01491 case kPiPlus:
01492 return 0.074085*std::pow(1.0-xf, 11.2413)*(1.0+168369.0*xf)*std::pow(xf, -1.0*2.47553);
01493 case kPiMinus:
01494 break;
01495 case kKPlus:
01496 break;
01497 case kKMinus:
01498 break;
01499 default:
01500 break;
01501 }
01502 }
01503 else
01504 {
01505 //
01506 //form: Par[4]*(1.0-xf)*std::exp(Par[5]*xf);
01507 //
01508 switch (particle)
01509 {
01510 case kPiPlus:
01511 return 49071.3*(1.0-xf)*std::exp(-12.1073*xf);
01512 case kPiMinus:
01513 break;
01514 case kKPlus:
01515 break;
01516 case kKMinus:
01517 break;
01518 default:
01519 break;
01520 }
01521 }
01522
01523
01524 return 0.0;
01525 }
|
|
||||||||||||
|
Definition at line 1528 of file Lfluk.cxx. References kKMinus, kKPlus, kPiMinus, and kPiPlus. Referenced by GetSKPCWeight(). 01529 {
01530 if (xf<0.1)
01531 {
01532 //
01533 //form: Par[0]*std::exp(Par[1]*std::pow(xf, Par[2]));
01534 //
01535 switch (particle)
01536 {
01537 case kPiPlus:
01538 return 35.8788*std::exp(-1.95477*std::pow(xf, 0.0638124));
01539 case kPiMinus:
01540 break;
01541 case kKPlus:
01542 break;
01543 case kKMinus:
01544 break;
01545 default:
01546 break;
01547 }
01548 }
01549 else
01550 {
01551 //
01552 //form: Par[3] + par[4]*xf + par[5]*xf*xf +par[6]*xf*xf*xf
01553 //
01554 switch (particle)
01555 {
01556 case kPiPlus:
01557 return 7.47623 + -7.35376*xf + -14.6949*xf*xf + 22.7689*xf*xf*xf;
01558 case kPiMinus:
01559 break;
01560 case kKPlus:
01561 break;
01562 case kKMinus:
01563 break;
01564 default:
01565 break;
01566 }
01567 }
01568
01569
01570 return 0.0;
01571
01572 }
|
|
||||||||||||
|
Definition at line 1575 of file Lfluk.cxx. References kKMinus, kKPlus, kPiMinus, and kPiPlus. Referenced by GetSKPCWeight(). 01576 {
01577 if (xf<0.25)
01578 {
01579 //
01580 //form: Par[0]*std::exp(Par[1]*std::pow(xf, Par[2]))*(Par{3] + Par[4]*xf + Par[5]*xf*xf);
01581 //
01582 switch (particle)
01583 {
01584 case kPiPlus:
01585 return 3.88469*std::exp(-6.20932*std::pow(xf, 0.26879))*(1.29961 + 53.0232*xf + 121.575*xf*xf);
01586 case kPiMinus:
01587 break;
01588 case kKPlus:
01589 break;
01590 case kKMinus:
01591 break;
01592 default:
01593 break;
01594 }
01595 }
01596 else
01597 {
01598 //
01599 //form: (Par[6]+Par[7]*xf)/std::pow(xf, Par[8]);
01600 //
01601 switch (particle)
01602 {
01603 case kPiPlus:
01604 return (-0.26414 + 1.87321*xf)/std::pow(xf, 1.28207);
01605 case kPiMinus:
01606 break;
01607 case kKPlus:
01608 break;
01609 case kKMinus:
01610 break;
01611 default:
01612 break;
01613 }
01614 }
01615
01616
01617 return 0.0;
01618
01619 }
|
|
||||||||||||||||
|
Definition at line 646 of file Lfluk.cxx. References C, GetSKPCA(), GetSKPCB(), GetSKPCC(), GetSKPCD(), GetSKPCE(), and SkewSKPC(). Referenced by CalcWeight(). 00647 {
00648 double weight = 1.0;
00649 double A, B, C, D, E;
00650 double Ap, Bp, Cp, Dp, Ep;
00651
00652 if (particle==kPiPlus)
00653 {
00654 if (pT>=0.03)
00655 {
00656 A = Lfluk::GetSKPCA(particle, xF);
00657 B = Lfluk::GetSKPCB(particle, xF);
00658 C = Lfluk::GetSKPCC(particle, xF);
00659 D = Lfluk::GetSKPCD(particle, xF);
00660 E = Lfluk::GetSKPCE(particle, xF);
00661
00662 Lfluk::SkewSKPC(xF, A, B, C, D, E, Ap, Bp, Cp, Dp, Ep);
00663
00664 if(Ep < -20.0) Ep = -20.0;
00665 double arg = -Dp*std::pow(pT,Ep)+D*std::pow(pT,E);
00666 if(arg > 100.0) arg = 100.0;
00667 if(arg < -100.0) arg = -100.0;
00668 weight=(Ap+Bp*pT+Cp*pT*pT)/(A+B*pT+C*pT*pT)*std::exp(arg);
00669 }
00670 else
00671 {
00672 weight = Lfluk::GetSKPCWeight(particle, 0.031, xF); // for low pT
00673 }
00674 }
00675 else if (particle==kKPlus)
00676 {
00677 weight = 1.0;
00678 }
00679 else if (particle==kPiMinus)
00680 {
00681 weight = 1.0;
00682 }
00683 else if (particle==kKMinus)
00684 {
00685 weight = 1.0;
00686 }
00687
00688 return weight;
00689 }
|
|
||||||||||||
|
Definition at line 1050 of file Lfluk.cxx. References kKMinus, kKPlus, kPiMinus, and kPiPlus. Referenced by GetSKZPWeight(). 01051 {
01052
01053 switch (particle)
01054 {
01055 case kPiPlus:
01056 return -0.007607*std::pow(1.0-xf, 4.045)*(1.0+0.9620e+4*xf)*std::pow(xf, -2.975);
01057 case kPiMinus:
01058 return -0.006306*std::pow(1.0-xf, 5.730)*(1.0+0.1365e+5*xf)*std::pow(xf, -2.900);
01059 case kKPlus:
01060 return -0.005187*std::pow(1.0-xf, 4.119)*(1.0+0.2170e+4*xf)*std::pow(xf, -2.767);
01061 case kKMinus:
01062 return -8.854 *std::pow(1.0-xf, 6.778)*(1.0-0.6050 *xf)*std::pow(xf, -1.827);
01063 default:
01064 break;
01065 }
01066
01067 return 1.0;
01068 }
|
|
||||||||||||
|
Definition at line 1071 of file Lfluk.cxx. References kKMinus, kKPlus, kPiMinus, and kPiPlus. Referenced by GetSKZPWeight(). 01072 {
01073 switch (particle)
01074 {
01075 case kPiPlus:
01076 return 0.05465*std::pow(1.0-xf, 2.675)*(1.0+0.6959e+5*xf)*std::pow(xf, -3.144);
01077 case kPiMinus:
01078 return 0.04608*std::pow(1.0-xf, 3.291)*(1.0+0.5857e+5*xf)*std::pow(xf, -3.209);
01079 case kKPlus:
01080 return 0.4918 *std::pow(1.0-xf, 2.672)*(1.0+0.1373e+4*xf)*std::pow(xf, -2.927);
01081 case kKMinus:
01082 return 0.02857*std::pow(1.0-xf, 7.494)*(1.0+0.5879e+5*xf)*std::pow(xf, -2.577);
01083 default:
01084 break;
01085 }
01086
01087 return 0.0;
01088 }
|
|
||||||||||||
|
Definition at line 1091 of file Lfluk.cxx. References kKMinus, kKPlus, kPiMinus, and kPiPlus. Referenced by GetSKZPWeight(). 01092 {
01093 if (xf<0.22)
01094 {
01095 switch (particle)
01096 {
01097 case kPiPlus:
01098 return -7.058/std::pow(xf, -0.1419 ) + 9.188;
01099 case kPiMinus:
01100 return -16.52 /std::pow(xf, -0.06204) + 18.12;
01101 case kKPlus:
01102 return -16.10 /std::pow(xf, -0.04582) + 17.92;
01103 case kKMinus:
01104 return -16.13 /std::pow(xf, -0.05678) + 17.39;
01105 default:
01106 break;
01107 }
01108 }
01109 else
01110 {
01111 switch (particle)
01112 {
01113 case kPiPlus:
01114 return 3.008/std::exp((xf-0.1948)*3.577) + 2.616 *xf + 0.1225;
01115 case kPiMinus:
01116 return 2.972/std::exp((xf-0.1758)*2.266) + 1.730 *xf + 0.04196;
01117 case kKPlus:
01118 return 6.905/std::exp((xf+0.1630)*6.718) - 0.4257*xf + 2.486;
01119 case kKMinus:
01120 return 3.916/std::exp((xf-4.615) *0.03255) - 4.702 *xf + 0.4062;
01121 default:
01122 break;
01123 }
01124 }
01125
01126
01127 return 0.0;
01128 }
|
|
||||||||||||||||
|
Definition at line 600 of file Lfluk.cxx. References C, GetSKZPA(), GetSKZPB(), GetSKZPC(), and SkewSKZP(). Referenced by CalcWeight(). 00601 {
00602 double weight = 1.0;
00603 double A,B,C;
00604 long double Ap,Bp,Cp;
00605
00606 if (particle==kPiPlus)
00607 {
00608 //Calculate weight for pions
00609 if (pT>=0.03)
00610 {
00611 // calculate the A, B and C in SKZP parameterization
00612 // A, B and C are best fit to Fluka 05
00613
00614 A = Lfluk::GetSKZPA(particle, xF);
00615 B = Lfluk::GetSKZPB(particle, xF);
00616 C = Lfluk::GetSKZPC(particle, xF);
00617
00618 Lfluk::SkewSKZP(xF, A, B, C, Ap, Bp, Cp);
00619
00620 if(Cp < -100.0) Cp = -100.0;
00621 if(Cp > 100.0) Cp = 100.0;
00622 weight=(Ap+Bp*pT)/(A+B*pT)*exp(-(Cp-C)*pow(pT,3./2.));
00623 }
00624 else
00625 {
00626 weight = Lfluk::GetSKZPWeight(particle, 0.031, xF); // for low pT
00627 }
00628 }
00629 else if (particle==kKPlus)
00630 {
00631 weight = 1.0;
00632 }
00633 else if (particle==kPiMinus)
00634 {
00635 weight = 1.0;
00636 }
00637 else if (particle==kKMinus)
00638 {
00639 weight = 1.0;
00640 }
00641
00642 return weight;
00643 }
|
|
||||||||||||||||
|
Definition at line 463 of file Lfluk.cxx. References CalcWeight(), fWeightHist, GetParticleEnum(), GetParticleName(), kKPlus, kPiMinus, kPiPlus, and ParticleType_t. 00464 {
00465
00466
00467 if(fWghtFromHist)
00468 {
00469
00470 double weight = 1.0;
00471
00472 Lfluk::ParticleType_t eptype=GetParticleEnum(ptype);
00473 string name = GetParticleName(eptype);
00474
00475 if(eptype != kPiPlus && eptype != kKPlus && eptype != kPiMinus && eptype != kPiMinus)
00476 {
00477 return 1.0;
00478 }
00479
00480 if (pz > 120.0 || pz < 0.0 || pT < 0.0)
00481 {
00482 return weight;
00483 }
00484
00485 const int nbinsx = fWeightHist[eptype] -> GetNbinsX();
00486 const int nbinsy = fWeightHist[eptype] -> GetNbinsY();
00487
00488 int binx = (fWeightHist[eptype] -> GetXaxis()) -> FindFixBin(pz);
00489 int biny = (fWeightHist[eptype] -> GetYaxis()) -> FindFixBin(pT);
00490
00491 if (pT > 1.0) biny = nbinsy+1;
00492
00493 if(!(binx < 0) && !(biny < 0) && !(binx > nbinsx+1) && !(biny > nbinsy+1))
00494 {
00495 weight = fWeightHist[eptype] -> GetBinContent(binx, biny);
00496 }
00497
00498 /*
00499 double wt = CalcWeight(eptype,pT,pz);
00500 string name = GetParticleName(eptype);
00501 //if(name == "PiMinus" && pz > 117 && pz < 119 && weight-wt != 0)
00502 if(weight - wt != 0)
00503 {
00504 cout << " Lfluk::GetWeight(): ptype = " << name << " pz = " << pz << " pt = " << pT << " weight from calc = "
00505 << setprecision(10) << wt << " weight from hist = " << weight << " diff = " << weight - wt << endl;
00506 }
00507 */
00508
00509 /*
00510 if((name == "PiPlus" || name == "KPlus") && (binx == 2) && biny == nbinsy+1)
00511 {
00512 cout << "get weight called : name = " << name << " pz = " << pz << " pt = " << pT << " weight = " << weight << endl;
00513 }
00514
00515 */
00516
00517
00518 return weight;
00519 }
00520
00521 else
00522 {
00523 return CalcWeight(ptype,pT,pz);
00524 }
00525 }
|
|
|
Definition at line 452 of file Lfluk.cxx. References fWeightHist, and GetParticleEnum(). 00453 {
00454 TH2* hist = 0;
00455 hist = dynamic_cast <TH2*> (fWeightHist[GetParticleEnum(ptype)]->Clone());
00456 hist->SetDirectory(0);
00457
00458 return hist;
00459
00460 }
|
|
|
Definition at line 45 of file Lfluk.h. 00045 {return fWeightHist[GetParticleEnum(ptype)];};
|
|
|
Definition at line 43 of file Lfluk.cxx. References base, fFluka05KPiRatio, fFluka05PiRatio, fFlukKPiRat, fFlukPiRat, fKPiRatio, fMeanPT, fN, fNA49Fit, fNA49pirat, fNWeighted, Form(), fPath, fPiRatio, fPlist, fPTXF, fWeightedMeanPT, fWeightedPTXF, fWeightHist, GetParticleName(), GetPTXF(), gSystem(), isFirst, kK0L, kKMinus, kKPlus, kPiMinus, kPiPlus, and MSG. 00044 {
00045 cout << "Lfluk::Init" << endl;
00046
00047 std::string dir = "";
00048 TH1::AddDirectory(kFALSE);
00049
00050 // topDir is where the fluka05 pt-xf histogram file is
00051 std::string topDir=dir; // user may set location of input data
00052 if(topDir=="") { // by default, this code looks in a standard place
00053 topDir="MCReweight/data";
00054 std::string base="";
00055 base=getenv("SRT_PRIVATE_CONTEXT");
00056 if (base!="" && base!=".") {
00057 // check if directory exists in SRT_PRIVATE_CONTEXT
00058 std::string path = base + "/" + topDir;
00059 void *dir_ptr = gSystem -> OpenDirectory(path.c_str());
00060 if(!dir_ptr) base=getenv("SRT_PUBLIC_CONTEXT"); // if it doesn't exist then use SRT_PUBLIC_CONTEXT
00061 }
00062 else base=getenv("SRT_PUBLIC_CONTEXT");
00063 if(base=="") {
00064 MSG("MCReweight",Msg::kFatal)<<"No SRT_PUBLIC_CONTEXT set"<<std::endl;
00065 assert(false);
00066 }
00067 topDir = base+ "/" + topDir;
00068 }
00069
00070 MSG("Lfluk",Msg::kDebug) <<"Lfluk reading data from: "<<topDir<<std::endl;
00071 std::string fileName =topDir+"/fluka05ptxf.root";
00072
00073 TFile* hFile=new TFile(fileName.c_str(), "READ");
00074 if(!(hFile->IsOpen()))
00075 {
00076 hFile = new TFile("/minos/scratch/loiacono/Data/ForMonFit/fluka05ptxf.root", "READ");
00077 if(!hFile || !(hFile->IsOpen()))
00078 {
00079 cout << "Can't open /minos/scratch/loiacono/Data/ForMonFit/fluka05ptxf.root" << endl;
00080 cout << " Failed to open file: " << fileName << endl;
00081 cout << " **************I am going to crash***************"<< endl;
00082 }
00083
00084 }
00085 else
00086 {
00087 cout << " Reading Fluka PTXF from: " << fileName << endl;
00088 }
00089
00090
00091 fPlist.push_back(kPiPlus);
00092 fPlist.push_back(kPiMinus);
00093 fPlist.push_back(kKPlus);
00094 fPlist.push_back(kKMinus);
00095 fPlist.push_back(kK0L);
00096
00097 for (std::vector<ParticleType_t>::iterator itPlist=fPlist.begin();itPlist!=fPlist.end(); itPlist++)
00098 {
00099 TH2F* hist=dynamic_cast <TH2F*> (hFile->Get(Form("hF05ptxf%s",GetParticleName(*itPlist).c_str()))->Clone());
00100 hist->SetDirectory(0);
00101 TH2F* hist2=dynamic_cast <TH2F*> (hist->Clone(Form("hWeightedPTXF%s",GetParticleName(*itPlist).c_str())));
00102 hist2->SetDirectory(0);
00103 hist2->SetTitle(Form("%s weighted pt-pz",GetParticleName(*itPlist).c_str()));
00104
00105 const int nbinsx = 120;
00106 const int nbinsy = hist -> GetNbinsY();
00107 //cout << "nbinsx = " << nbinsx << "nbinsy = " << nbinsy << endl;
00108 //fWeightHist[*itPlist]=dynamic_cast<TH2D*> (hist->Clone(Form("hWeight%s",GetParticleName(*itPlist).c_str())));
00109 //
00110 //The bins need to be the same as the bins that the pt-pz plots for the horn current scans!!!!!!!!!!!!!!!!!!!!!!!
00111 //
00112 fWeightHist[*itPlist] = new TH2D(Form("hWeightedPTXF%s",GetParticleName(*itPlist).c_str()),GetParticleName(*itPlist).c_str(), nbinsx, 0, 120, nbinsy, 0, 1);
00113 for(int ibinx = 0; ibinx <= nbinsx+1; ++ibinx)
00114 {
00115 for(int ibiny = 0; ibiny <= nbinsy+1; ++ibiny)
00116 {
00117 //cout << "ibinx = " << ibinx << "ibiny = " << ibiny << endl;
00118 fWeightHist[*itPlist] -> SetBinContent(ibinx, ibiny, 1.0);
00119 }
00120 }
00121 //TH2D* temp = dynamic_cast<TH2D*>(hist->Clone());
00122 //fWeightHist[*itPlist]->Divide(temp);
00123 fWeightHist[*itPlist]->SetDirectory(0);
00124
00125 std::pair<Lfluk::ParticleType_t, TH2F* > p(*itPlist, hist);
00126 std::pair<Lfluk::ParticleType_t, TH2F* > p2(*itPlist, hist2);
00127 fPTXF.insert(p);
00128 fWeightedPTXF.insert(p2);
00129
00130 fPTXF[*itPlist]->SetDirectory(0);
00131 fWeightedPTXF[*itPlist]->SetDirectory(0);
00132
00133 //
00134 // Projection DOES include overflow bins
00135 // Get Mean DOES NOT
00136 //
00137 fMeanPT[*itPlist]=fPTXF[*itPlist]->ProjectionY()->GetMean()*1000.;
00138 fWeightedMeanPT[*itPlist]=fWeightedPTXF[*itPlist]->ProjectionY()->GetMean()*1000.;
00139
00140 double N=fPTXF[*itPlist]->ProjectionY()->GetSumOfWeights();
00141 double wN=fWeightedPTXF[*itPlist]->ProjectionY()->GetSumOfWeights();
00142
00143 fN[*itPlist]=N;
00144 fNWeighted[*itPlist]=wN;
00145 }
00146 hFile->Close();
00147 delete hFile;
00148 hFile=0;
00149
00150 //
00151 //Get NA49 Ratio from Fit.
00152 //
00153 std::string path = string(fPath + "na49Ratio.root");
00154 TFile *fitfile = new TFile(path.c_str(), "READ");
00155 if(!(fitfile->IsOpen()))
00156 {
00157 cout << " Failed to open file: " << path << endl;
00158 cout << " **************I am going to crash***************"<< endl;
00159 }
00160 else
00161 {
00162 cout << " Reading NA49 Pi Ratio from: " << path << endl;
00163 }
00164
00165 stringstream name;
00166 name << "NA49PipPimRatio_" << fNA49Fit << "_fit";
00167 if(!(fitfile -> Get(name.str().c_str())))
00168 {
00169 cout << " Failed to get " << name.str() << " from " << path << endl;
00170 }
00171 fNA49pirat = dynamic_cast<TF1 *>( (fitfile -> Get(name.str().c_str())) -> Clone());
00172 if(!fNA49pirat)
00173 {
00174 cout << " Failed to get " << name.str() << " from " << path << endl;
00175 fNA49pirat = 0;
00176 }
00177 fitfile -> Close();
00178 delete fitfile;
00179
00180
00181 /*
00182 //
00183 //added stuff above xf=0.5.
00184 //see MuonPhysics/Fit/Misc/NA49PiRatio.cxx
00185 //
00186
00187 //na49 data thru xf = 0.5
00188 double xfna49[] ={0.01, 0.02, 0.03, 0.04, 0.05,
00189 0.075, 0.1, 0.125, 0.15, 0.2,
00190 0.25, 0.3, 0.4, 0.5, 70.5,
00191 80.5, 90.5, 100.5, 119.5};
00192 double pipna49[]={25.756, 20.903, 16.683, 13.389, 10.981,
00193 7.067, 4.904, 3.531, 2.648, 1.5499,
00194 0.9734, 0.5996, 0.2455, 0.1081, 3.69741,
00195 4.36113,5.3943, 5.89819,7.60617};
00196 double pimna49[]={24.424, 19.247, 14.938, 11.757, 9.516,
00197 5.790, 3.803, 2.561, 1.826, 0.9555,
00198 0.5475, 0.3190, 0.1041, 0.0336, 1.0,
00199 1.0, 1.0, 1.0, 1.0};
00200 //double piratna49[14];
00201 double piratna49[19];
00202
00203 for (int i=0;i<19;i++) piratna49[i]=pipna49[i]/pimna49[i];
00204 for (int i=0;i<14;i++) xfna49[i]*=120.;
00205 TGraph* pina49=new TGraph(19,xfna49,piratna49);
00206 fNA49pirat = new TMVA::TSpline1("sp0", pina49);
00207 */
00208
00209 fPiRatio.resize(122);
00210 fFluka05KPiRatio.resize(122);
00211 fFluka05PiRatio.resize(122);
00212 fKPiRatio.resize(122);
00213 TH1D* piplus=dynamic_cast<TH1D*> (GetPTXF(8)
00214 ->ProjectionX()->Clone("piplus"));
00215 TH1D* piminus=dynamic_cast<TH1D*> (GetPTXF(9)
00216 ->ProjectionX()->Clone("piplus"));
00217 TH1D* kplus=dynamic_cast<TH1D*> (GetPTXF(11)
00218 ->ProjectionX()->Clone("kplus"));
00219
00220 kplus->Divide(piplus);
00221 piplus->Divide(piminus);
00222
00223 fFlukKPiRat = dynamic_cast<TH1D*>( kplus->Clone("FlukKPiRatio") );
00224 fFlukPiRat = dynamic_cast<TH1D*>( piplus->Clone("FlukPiRatio") );
00225
00226
00227
00228 for (int i=1;i<121;i++)
00229 {
00230 if (kplus->GetBinContent(i)!=0.)
00231 fFluka05KPiRatio[i] = kplus->GetBinContent(i);
00232 else if (i<119)
00233 fFluka05KPiRatio[i] = (fFluka05KPiRatio[i-2]+fFluka05KPiRatio[i-1]
00234 +kplus->GetBinContent(i+1)+kplus->GetBinContent(i+2))/4.;
00235 else
00236 fFluka05KPiRatio[i] = fFluka05KPiRatio[i-1];
00237
00238 fFlukKPiRat -> SetBinContent(i, fFluka05KPiRatio[i]);
00239 }
00240
00241
00242
00243 for (int i=1;i<121;i++)
00244 {
00245 if (piplus->GetBinContent(i)!=0.)
00246 fFluka05PiRatio[i] = piplus->GetBinContent(i);
00247 else if (i<119)
00248 fFluka05PiRatio[i] = (fFluka05PiRatio[i-2]+fFluka05PiRatio[i-1]
00249 +piplus->GetBinContent(i+1)+piplus->GetBinContent(i+2))/4.;
00250 else
00251 fFluka05PiRatio[i] = fFluka05PiRatio[i-1];
00252
00253 fFlukPiRat -> SetBinContent(i, fFluka05PiRatio[i]);
00254 }
00255
00256
00257 isFirst=true;
00258
00259
00260 }
|
|
|
Definition at line 1132 of file Lfluk.cxx. References CalcWeight(), fFluka05KPiRatio, fFluka05PiRatio, fKPiRatio, fNA49pirat, fNWeighted, fPiRatio, fPlist, fPTXF, fUseFlukPi, fWeightedMeanPT, fWeightedPTXF, fWeightHist, GetParticleName(), GetReweightedPTXF(), isFirst, and ParticleType_t. Referenced by SetParameters(). 01133 {
01134
01135 if (isFirst) {
01136 for (int i=0;i<121;i++) {
01137 fKPiRatio[i]=1.;
01138 fPiRatio[i]=1.;
01139 }
01140 }
01141 for (std::vector<ParticleType_t>::iterator itPlist=fPlist.begin();
01142 itPlist!=fPlist.end(); itPlist++)
01143 {
01144 double xf,pt;
01145 double meanpt(0.);
01146 double sum(0.);
01147 const int nbinsx = fPTXF[*itPlist]->GetNbinsX();
01148 const int nbinsy = fPTXF[*itPlist]->GetNbinsY();
01149 for (int i=0;i<=nbinsx+1;i++)
01150 {
01151 for (int j=0;j<=nbinsy+1;j++)
01152 {
01153 xf=fPTXF[*itPlist]->GetXaxis()->GetBinCenter(i);
01154 pt=fPTXF[*itPlist]->GetYaxis()->GetBinCenter(j);
01155 double weight = CalcWeight(*itPlist,pt,xf);
01156
01157 fWeightedPTXF[*itPlist]
01158 ->SetBinContent(i,j,(fPTXF[*itPlist]->GetBinContent(i,j)
01159 *weight));
01160
01161 Lfluk::ParticleType_t eptype= *itPlist;
01162 string name = GetParticleName(eptype);
01163 fWeightHist[*itPlist]->SetBinContent(i,j,weight);
01164 if(weight - fWeightHist[*itPlist]->GetBinContent(i,j) != 0)
01165 //if(name == "PiPlus" && weight !=1.0)
01166 {
01167 cout << "Lfluk::RecalculateWeights(): name = " << name << " pz = " << xf << " pt = " << pt
01168 << " weight from calc = " << setprecision(10) << weight << " weight from hist = "
01169 << fWeightHist[*itPlist]->GetBinContent(i,j)
01170 << " diff = " << weight - fWeightHist[*itPlist]->GetBinContent(i,j) << endl;
01171 }
01172
01173 //
01174 // Projection DOES include overflow bins
01175 // Get Mean DOES NOT
01176 //
01177 //need to include the overflow in y for sure in fWeightedPTXF
01178 //because it is used in the ratios below
01179 //but in meanpt we don't want to use this bin
01180 //I have checked that the fluka "signal" in i=nbinsx+1 is always 0
01181 //so it doesn't matter if this bin is included or not. I think the
01182 //proper thing is that is should be because the projection does
01183 //include it
01184 //
01185 if(j == nbinsy+1) continue;
01186 meanpt+=fWeightedPTXF[*itPlist]->GetBinContent(i,j)*pt;
01187 sum+=fWeightedPTXF[*itPlist]->GetBinContent(i,j);
01188
01189 /*Lfluk::ParticleType_t eptype= *itPlist;
01190 string name = GetParticleName(eptype);
01191 if((name == "PiPlus" || name == "KPlus") && ((i == 2 && j == nbinsy+1) || (i==nbinsx+1 && j == 2)))
01192 {
01193 cout <<"name = " << name << " pz = " << xf << " pt = " << pt
01194 << " weight = " << weight << " fluksig = " << fPTXF[*itPlist]->GetBinContent(i,j)
01195 << " newsig = " << fWeightedPTXF[*itPlist]->GetBinContent(i,j) << endl;
01196 }*/
01197
01198
01199 }
01200 }
01201
01202 meanpt/=sum;
01203 fWeightedMeanPT[*itPlist]=meanpt*1000.; //GeV to MeV
01204 fNWeighted[*itPlist]=sum;
01205 meanpt=0.;
01206 sum=0.;
01207 /*
01208 if(!isFirst)
01209 {
01210 for (int i=0;i<=fWeightHist[*itPlist]->GetNbinsX()+1;i++)
01211 {
01212 for (int j=0;j<=fWeightHist[*itPlist]->GetNbinsY()+1;j++)
01213 {
01214 xf=fWeightHist[*itPlist]->GetXaxis()->GetBinCenter(i);
01215 pt=fWeightHist[*itPlist]->GetYaxis()->GetBinCenter(j);
01216 double weight = CalcWeight(*itPlist,pt,xf);
01217
01218 Lfluk::ParticleType_t eptype= *itPlist;
01219 string name = GetParticleName(eptype);
01220 fWeightHist[*itPlist]->SetBinContent(i,j,weight);
01221
01222
01223 //if(name == "PiMinus" && weight - fWeightHist[*itPlist]->GetBinContent(i,j) != 0)
01224 //if(name == "PiPlus" && xf > 49 && xf < 50 && weight - fWeightHist[*itPlist]->GetBinContent(i,j) != 0)
01225 //if(name == "PiMinus" && xf > 117 && xf < 119 && weight - fWeightHist[*itPlist]->GetBinContent(i,j) != 0)
01226 if(weight - fWeightHist[*itPlist]->GetBinContent(i,j) != 0)
01227 {
01228 cout << "Lfluk::RecalculateWeights(): name = " << name << " pz = " << xf << " pt = " << pt
01229 << " weight from calc = " << setprecision(10) << weight << " weight from hist = "
01230 << fWeightHist[*itPlist]->GetBinContent(i,j)
01231 << " diff = " << weight - fWeightHist[*itPlist]->GetBinContent(i,j) << endl;
01232 }
01233 }
01234 }
01235 }
01236 */
01237 }
01238
01239
01240 if (isFirst)
01241 {
01242 isFirst=false;
01243 const TH1D* piplus=dynamic_cast<TH1D*> (GetReweightedPTXF(8)
01244 ->ProjectionX()->Clone("piplus"));
01245 const TH1D* piminus=dynamic_cast<TH1D*> (GetReweightedPTXF(9)
01246 ->ProjectionX()->Clone("piminus"));
01247 TH1D* kplus_f05=dynamic_cast<TH1D*> (GetReweightedPTXF(11)
01248 ->ProjectionX()->Clone("kplus_f05"));
01249 TH1D* piplus_rat=dynamic_cast<TH1D*> (GetReweightedPTXF(8)
01250 ->ProjectionX()->Clone("piplus_f05"));
01251
01252 //fPiPlusReWProj = dynamic_cast<TH1D*>(piplus->Clone());
01253
01254 piplus_rat -> Divide(piplus, piminus);
01255
01256 //fKPlusf05Proj = dynamic_cast<TH1D*>(kplus_f05->Clone());
01257
01258 kplus_f05->Divide(piplus);
01259
01260 //fKPlusf05PipReWRatio = dynamic_cast<TH1D*>(kplus_f05->Clone());
01261
01262 for (int i=1;i<121;i++)
01263 {
01264 if (kplus_f05->GetBinContent(i)!=0)
01265 {
01266 fKPiRatio[i] = (fFluka05KPiRatio[i]*fKPiRatScale)/kplus_f05->GetBinContent(i);
01267 }
01268
01269 if (piplus_rat->GetBinContent(i)!=0)
01270 {
01271 if(fUseFlukPi && fFluka05PiRatio[i] != 0) fPiRatio[i] = piplus_rat->GetBinContent(i)/fFluka05PiRatio[i];
01272 else fPiRatio[i] = piplus_rat->GetBinContent(i)/fNA49pirat->Eval(piplus_rat->GetBinCenter(i));
01273
01274 }
01275 /*if (i<62)
01276 {
01277 fPiRatio[i] = piplus->GetBinContent(i)/fNA49pirat->Eval(piplus->GetBinCenter(i));
01278 }
01279 else
01280 {
01281 fPiRatio[i]=1.;
01282 }*/
01283 }
01284 delete piplus;
01285 piplus=0;
01286 delete piminus;
01287 piminus=0;
01288 delete kplus_f05;
01289 kplus_f05=0;
01290 delete piplus_rat;
01291 piplus_rat = 0;
01292
01293 RecalculateWeights();
01294 }
01295
01296 }
|
|
||||||||||||||||||||||||
|
Definition at line 262 of file Lfluk.cxx. References fMaxW, fMinW, fNA49Fit, fParam, and fSkew. 00263 {
00264 fParam = param;
00265 fSkew = skew;
00266 fNA49Fit = na49fit;
00267 fMinW = wmin;
00268 fMaxW = wmax;
00269
00270 }
|
|
|
Definition at line 352 of file Lfluk.cxx. References fPar, fParam, fSkew, isFirst, and RecalculateWeights(). 00353 {
00354 fPar=par;
00355 isFirst=true;
00356
00357 if(fParam == "skzp" && fSkew == "fluka")
00358 {
00359 cout << " Pars = ";
00360 for(unsigned int ip = 0; ip < fPar.size(); ++ip)
00361 {
00362 cout << fPar[ip] << " ";
00363 }
00364 cout << "end pars " << endl;
00365
00366 /*
00367 Lfluk::ParticleType_t eptype=GetParticleEnum(8);
00368
00369 double xF; bool print = false;
00370 for (int i=0;i<fPTXF[eptype]->GetNbinsX()+1;i++)
00371 {
00372 xF=fPTXF[eptype]->GetXaxis()->GetBinCenter(i);
00373 xF = xF/120.0;
00374 if (xF>1.||xF<0.) continue;
00375
00376 if( fabs((int)(xF*10) - xF*10) < 0.1) print = true;
00377
00378 if(print) cout << "xF = " << xF;
00379
00380 double A = Lfluk::GetSKZPA(eptype, xF);
00381 if(print) cout << " A = " << A;
00382
00383 double B = Lfluk::GetSKZPB(eptype, xF);
00384 if(print) cout << " B = " << B;
00385
00386 double C = Lfluk::GetSKZPC(eptype, xF);
00387 if(print) cout << " C = " << C;
00388
00389 long double Ap= fPar[0]*pow((1.-xF),fPar[1])*(1.+fPar[2]*xF)*pow(xF,fPar[3]);
00390 if(print) cout << " Ap = " << Ap;
00391
00392 long double Bp= fPar[4]*pow((1.-xF),fPar[5])*(1.+fPar[6]*xF)*pow(xF,fPar[7]);
00393 if(print) cout << " Bp = " << Bp;
00394
00395 long double Cp;
00396 if (xF<0.22)
00397 {
00398 Cp= fPar[8]/pow(xF,fPar[9])+fPar[10];
00399 if(Cp < -100.0) Cp = -100.0;
00400 if(Cp > 100.0) Cp = 100.0;
00401 if(print) cout << " Cp = " << Cp;
00402 }
00403 else
00404 {
00405 double arg22 = (0.22-fPar[11])*fPar[12];
00406 if(arg22 < -100.0) arg22 = -100.0;
00407 if(arg22 > 100.0) arg22 = 100.0;
00408 double match = (fPar[8]/pow(0.22,fPar[9])+fPar[10]-
00409 fPar[13]*0.22-fPar[14])*exp(arg22);
00410 double exparg = (xF-fPar[11])*fPar[12];
00411 if(exparg < -100.0) exparg = -100.0;
00412 if(exparg > 100.0) exparg = 100.0;
00413 Cp = ( match / exp(exparg)) + fPar[13]*xF+fPar[14];
00414
00415 if(Cp < -100.0) Cp = -100.0;
00416 if(Cp > 100.0) Cp = 100.0;
00417 if(print) cout << " Cp = " << Cp;
00418 }
00419
00420
00421 if(print)
00422 {
00423 cout << " exp(-(Cp-C)) = " << exp(-(Cp-C)) << endl;
00424
00425 cout << " pT = 0.01, (Ap+Bp*pT)/(A+B*pT) = " << (Ap+Bp*0.01)/(A+B*0.01);
00426 double weight=(Ap+Bp*0.01)/(A+B*0.01)*exp(-(Cp-C)*pow(0.01,3./2.));
00427 cout << " weight = " << weight << endl;
00428
00429 cout << " pT = 0.95, (Ap+Bp*pT)/(A+B*pT) = " << (Ap+Bp*0.95)/(A+B*0.95);
00430 weight=(Ap+Bp*0.95)/(A+B*0.95)*exp(-(Cp-C)*pow(0.95,3./2.));
00431 cout << " weight = " << weight << endl;
00432
00433 }
00434 print = false;
00435
00436 }
00437 */
00438
00439 cout << "Done Setting pars...recalc weights" << endl;
00440 cout << endl;
00441
00442 }
00443
00444
00445
00446
00447
00448 RecalculateWeights();
00449 }
|
|
||||||||||||||||||||||||||||||||||||||||||||||||
|
Definition at line 935 of file Lfluk.cxx. Referenced by GetSKPCWeight(). 00938 {
00939
00940 if(fSkew.find("lin") != string::npos)
00941 {
00942 if(fPar.size() == 10)
00943 {
00944 Ap=(fPar[0]+fPar[1]*xF)*A;
00945 Bp=(fPar[2]+fPar[3]*xF)*B;
00946 Cp=(fPar[4]+fPar[5]*xF)*C;
00947 Dp=(fPar[6]+fPar[7]*xF)*D;
00948 Ep=(fPar[8]+fPar[9]*xF)*E;
00949
00950 return;
00951 }
00952 }
00953
00954 if(fSkew.find("quad") != string::npos)
00955 {
00956 if(fPar.size() == 15)
00957 {
00958 Ap=(fPar[0]+fPar[1]*xF+fPar[2]*xF*xF)*A;
00959 Bp=(fPar[3]+fPar[4]*xF+fPar[5]*xF*xF)*B;
00960 Cp=(fPar[6]+fPar[7]*xF+fPar[8]*xF*xF)*C;
00961 Dp=(fPar[9]+fPar[10]*xF+fPar[11]*xF*xF)*D;
00962 Ep=(fPar[12]+fPar[13]*xF+fPar[14]*xF*xF)*E;
00963
00964 return;
00965 }
00966 }
00967
00968 cout << "Lfluk::SkewSKPC() - No valid skew option or wrong number of parameters" << endl
00969 << " fSkew = " << fSkew << endl
00970 << " NPars = " << fPar.size() << endl;
00971
00972 Ap = A;
00973 Bp = B;
00974 Cp = C;
00975 Dp = D;
00976 Ep = E;
00977
00978 }
|
|
||||||||||||||||||||||||||||||||
|
Definition at line 981 of file Lfluk.cxx. Referenced by GetSKZPWeight(). 00982 {
00983 // scale/skew A, B and C
00984
00985 if(fSkew.find("fluka") != string::npos)
00986 {
00987 if(fPar.size() == 15)
00988 {
00989 Ap= fPar[0]*pow((1.-xF),fPar[1])*(1.+fPar[2]*xF)*pow(xF,fPar[3]);
00990 Bp= fPar[4]*pow((1.-xF),fPar[5])*(1.+fPar[6]*xF)*pow(xF,fPar[7]);
00991
00992 if (xF<0.22)
00993 {
00994 Cp= fPar[8]/pow(xF,fPar[9])+fPar[10];
00995 }
00996 else
00997 {
00998 double arg22 = (0.22-fPar[11])*fPar[12];
00999 if(arg22 < -100.0) arg22 = -100.0;
01000 if(arg22 > 100.0) arg22 = 100.0;
01001 double match = (fPar[8]/pow(0.22,fPar[9])+fPar[10]-
01002 fPar[13]*0.22-fPar[14])*exp(arg22);
01003 double exparg = (xF-fPar[11])*fPar[12];
01004 if(exparg < -100.0) exparg = -100.0;
01005 if(exparg > 100.0) exparg = 100.0;
01006 Cp = ( match / exp(exparg)) + fPar[13]*xF+fPar[14];
01007
01008 }
01009
01010 return;
01011 }
01012 }
01013
01014 if(fSkew.find("lin") != string::npos)
01015 {
01016 if(fPar.size() == 6)
01017 {
01018 Ap=(fPar[0]+fPar[1]*xF)*A;
01019 Bp=(fPar[2]+fPar[3]*xF)*B;
01020 Cp=(fPar[4]+fPar[5]*xF)*C;
01021
01022 return;
01023 }
01024 }
01025
01026 if(fSkew.find("quad") != string::npos)
01027 {
01028 if(fPar.size() == 9)
01029 {
01030 Ap=(fPar[0]+fPar[1]*xF+fPar[2]*xF*xF)*A;
01031 Bp=(fPar[3]+fPar[4]*xF+fPar[5]*xF*xF)*B;
01032 Cp=(fPar[6]+fPar[7]*xF+fPar[8]*xF*xF)*C;
01033
01034 return;
01035 }
01036 }
01037
01038 cout << "Lfluk::SkewSKZP() - No valid skew option or wrong number of parameters" << endl
01039 << " fSkew = " << fSkew << endl
01040 << " NPars = " << fPar.size() << endl;
01041
01042 Ap = A;
01043 Bp = B;
01044 Cp = C;
01045
01046 }
|
|
|
Definition at line 125 of file Lfluk.h. Referenced by Init(), and RecalculateWeights(). |
|
|
Definition at line 126 of file Lfluk.h. Referenced by Init(), and RecalculateWeights(). |
|
|
Definition at line 131 of file Lfluk.h. Referenced by Init(). |
|
|
Definition at line 132 of file Lfluk.h. Referenced by Init(). |
|
|
Definition at line 127 of file Lfluk.h. Referenced by CalcWeight(), Init(), and RecalculateWeights(). |
|
|
Definition at line 105 of file Lfluk.h. Referenced by Config(). |
|
|
|
|
|
|
|
|
Definition at line 104 of file Lfluk.h. Referenced by Config(), and SetConfig(). |
|
|
Definition at line 116 of file Lfluk.h. Referenced by GetPTshift(), and Init(). |
|
|
Definition at line 103 of file Lfluk.h. Referenced by Config(), and SetConfig(). |
|
|
Definition at line 117 of file Lfluk.h. Referenced by CalcWeight(), and Init(). |
|
|
Definition at line 100 of file Lfluk.h. Referenced by Config(), Init(), and SetConfig(). |
|
|
Definition at line 123 of file Lfluk.h. Referenced by Init(), and RecalculateWeights(). |
|
|
|
|
|
|
|
|
Definition at line 118 of file Lfluk.h. Referenced by CalcWeight(), Init(), and RecalculateWeights(). |
|
|
Definition at line 110 of file Lfluk.h. Referenced by CalcWeight(), GetMeanPTWeight(), GetQuadWeight(), SetParameters(), SkewSKPC(), and SkewSKZP(). |
|
|
Definition at line 98 of file Lfluk.h. Referenced by CalcWeight(), Config(), GetQuadWeight(), SetConfig(), and SetParameters(). |
|
|
|
|
|
|
|
|
Definition at line 124 of file Lfluk.h. Referenced by CalcWeight(), Init(), and RecalculateWeights(). |
|
|
Definition at line 111 of file Lfluk.h. Referenced by Init(), and RecalculateWeights(). |
|
|
Definition at line 113 of file Lfluk.h. Referenced by GetNFrac(), GetPTshift(), Init(), and RecalculateWeights(). |
|
|
Definition at line 99 of file Lfluk.h. Referenced by Config(), SetConfig(), SetParameters(), SkewSKPC(), and SkewSKZP(). |
|
|
Definition at line 108 of file Lfluk.h. Referenced by Config(), and RecalculateWeights(). |
|
|
Definition at line 119 of file Lfluk.h. Referenced by GetPTshift(), Init(), and RecalculateWeights(). |
|
|
Definition at line 114 of file Lfluk.h. Referenced by GetNFrac(), GetPTshift(), Init(), and RecalculateWeights(). |
|
|
Definition at line 115 of file Lfluk.h. Referenced by GetWeight(), GetWeightHist(), Init(), and RecalculateWeights(). |
|
|
Definition at line 107 of file Lfluk.h. Referenced by Config(). |
|
|
Definition at line 129 of file Lfluk.h. Referenced by Init(), RecalculateWeights(), and SetParameters(). |
1.3.9.1