rivet is hosted by Hepforge, IPPP Durham
SLD_1999_S3743934.cc
Go to the documentation of this file.
00001 // -*- C++ -*-
00002 #include "Rivet/Analysis.hh"
00003 #include "Rivet/Projections/Beam.hh"
00004 #include "Rivet/Projections/FinalState.hh"
00005 #include "Rivet/Projections/UnstableFinalState.hh"
00006 #include "Rivet/Projections/ChargedFinalState.hh"
00007 #include "Rivet/Projections/InitialQuarks.hh"
00008 #include "Rivet/Projections/Thrust.hh"
00009 
00010 namespace Rivet {
00011 
00012 
00013   /// @brief SLD flavour-dependent fragmentation paper
00014   /// @author Peter Richardson
00015   class SLD_1999_S3743934 : public Analysis {
00016   public:
00017 
00018     /// Constructor
00019     SLD_1999_S3743934()
00020       : Analysis("SLD_1999_S3743934"),
00021         _SumOfudsWeights(0.), _SumOfcWeights(0.),
00022         _SumOfbWeights(0.),
00023         _multPiPlus(4,0.),_multKPlus(4,0.),_multK0(4,0.),
00024         _multKStar0(4,0.),_multPhi(4,0.),
00025         _multProton(4,0.),_multLambda(4,0.)
00026     {    }
00027 
00028 
00029     /// @name Analysis methods
00030     //@{
00031 
00032     void analyze(const Event& e) {
00033       // First, veto on leptonic events by requiring at least 4 charged FS particles
00034       const FinalState& fs = applyProjection<FinalState>(e, "FS");
00035       const size_t numParticles = fs.particles().size();
00036 
00037       // Even if we only generate hadronic events, we still need a cut on numCharged >= 2.
00038       if (numParticles < 2) {
00039         MSG_DEBUG("Failed ncharged cut");
00040         vetoEvent;
00041       }
00042       MSG_DEBUG("Passed ncharged cut");
00043       // Get event weight for histo filling
00044       const double weight = e.weight();
00045 
00046       // Get beams and average beam momentum
00047       const ParticlePair& beams = applyProjection<Beam>(e, "Beams").beams();
00048       const double meanBeamMom = ( beams.first.momentum().vector3().mod() +
00049                                    beams.second.momentum().vector3().mod() ) / 2.0;
00050       MSG_DEBUG("Avg beam momentum = " << meanBeamMom);
00051       int flavour = 0;
00052       const InitialQuarks& iqf = applyProjection<InitialQuarks>(e, "IQF");
00053 
00054       // If we only have two quarks (qqbar), just take the flavour.
00055       // If we have more than two quarks, look for the highest energetic q-qbar pair.
00056       /// @todo Can we make this based on hadron flavour instead?
00057       Particles quarks;
00058       if (iqf.particles().size() == 2) {
00059         flavour = abs( iqf.particles().front().pdgId() );
00060         quarks = iqf.particles();
00061       } else {
00062         map<int, Particle > quarkmap;
00063         foreach (const Particle& p, iqf.particles()) {
00064           if (quarkmap.find(p.pdgId()) == quarkmap.end()) quarkmap[p.pdgId()] = p;
00065           else if (quarkmap[p.pdgId()].momentum().E() < p.momentum().E()) quarkmap[p.pdgId()] = p;
00066         }
00067         double maxenergy = 0.;
00068         for (int i = 1; i <= 5; ++i) {
00069           double energy(0.);
00070           if (quarkmap.find( i) != quarkmap.end())
00071             energy += quarkmap[ i].momentum().E();
00072           if (quarkmap.find(-i) != quarkmap.end())
00073             energy += quarkmap[-i].momentum().E();
00074           if (energy > maxenergy)
00075             flavour = i;
00076         }
00077         if (quarkmap.find(flavour) != quarkmap.end())
00078           quarks.push_back(quarkmap[flavour]);
00079         if (quarkmap.find(-flavour) != quarkmap.end())
00080           quarks.push_back(quarkmap[-flavour]);
00081       }
00082       switch (flavour) {
00083       case PID::DQUARK:
00084       case PID::UQUARK:
00085       case PID::SQUARK:
00086         _SumOfudsWeights += weight;
00087         break;
00088       case PID::CQUARK:
00089         _SumOfcWeights += weight;
00090         break;
00091       case PID::BQUARK:
00092         _SumOfbWeights += weight;
00093         break;
00094       }
00095       // thrust axis for projections
00096       Vector3 axis = applyProjection<Thrust>(e, "Thrust").thrustAxis();
00097       double dot(0.);
00098       if (!quarks.empty()) {
00099         dot = quarks[0].momentum().vector3().dot(axis);
00100         if (quarks[0].pdgId() < 0) dot *= -1;
00101       }
00102 
00103       foreach (const Particle& p, fs.particles()) {
00104         const double xp = p.momentum().vector3().mod()/meanBeamMom;
00105         // if in quark or antiquark hemisphere
00106         bool quark = p.momentum().vector3().dot(axis)*dot > 0.;
00107         _h_XpChargedN->fill(xp, weight);
00108         _temp_XpChargedN1->fill(xp, weight);
00109         _temp_XpChargedN2->fill(xp, weight);
00110         _temp_XpChargedN3->fill(xp, weight);
00111         int id = abs(p.pdgId());
00112         // charged pions
00113         if (id == PID::PIPLUS) {
00114           _h_XpPiPlusN->fill(xp, weight);
00115           _multPiPlus[0] += weight;
00116           switch (flavour) {
00117           case PID::DQUARK:
00118           case PID::UQUARK:
00119           case PID::SQUARK:
00120             _multPiPlus[1] += weight;
00121             _h_XpPiPlusLight->fill(xp, weight);
00122             if( ( quark && p.pdgId()>0 ) || ( !quark && p.pdgId()<0 ))
00123               _h_RPiPlus->fill(xp, weight);
00124             else
00125               _h_RPiMinus->fill(xp, weight);
00126             break;
00127           case PID::CQUARK:
00128             _multPiPlus[2] += weight;
00129             _h_XpPiPlusCharm->fill(xp, weight);
00130             break;
00131           case PID::BQUARK:
00132             _multPiPlus[3] += weight;
00133             _h_XpPiPlusBottom->fill(xp, weight);
00134             break;
00135           }
00136         }
00137         else if (id == PID::KPLUS) {
00138           _h_XpKPlusN->fill(xp, weight);
00139           _multKPlus[0] += weight;
00140           switch (flavour) {
00141           case PID::DQUARK:
00142           case PID::UQUARK:
00143           case PID::SQUARK:
00144             _multKPlus[1] += weight;
00145             _temp_XpKPlusLight->fill(xp, weight);
00146             _h_XpKPlusLight->fill(xp, weight);
00147             if( ( quark && p.pdgId()>0 ) || ( !quark && p.pdgId()<0 ))
00148               _h_RKPlus->fill(xp, weight);
00149             else
00150               _h_RKMinus->fill(xp, weight);
00151             break;
00152          break;
00153           case PID::CQUARK:
00154             _multKPlus[2] += weight;
00155             _h_XpKPlusCharm->fill(xp, weight);
00156             _temp_XpKPlusCharm->fill(xp, weight);
00157             break;
00158           case PID::BQUARK:
00159             _multKPlus[3] += weight;
00160             _h_XpKPlusBottom->fill(xp, weight);
00161             break;
00162           }
00163         }
00164         else if (id == PID::PROTON) {
00165           _h_XpProtonN->fill(xp, weight);
00166           _multProton[0] += weight;
00167           switch (flavour) {
00168           case PID::DQUARK:
00169           case PID::UQUARK:
00170           case PID::SQUARK:
00171             _multProton[1] += weight;
00172             _temp_XpProtonLight->fill(xp, weight);
00173             _h_XpProtonLight->fill(xp, weight);
00174             if( ( quark && p.pdgId()>0 ) || ( !quark && p.pdgId()<0 ))
00175               _h_RProton->fill(xp, weight);
00176             else
00177               _h_RPBar  ->fill(xp, weight);
00178             break;
00179          break;
00180           case PID::CQUARK:
00181             _multProton[2] += weight;
00182             _temp_XpProtonCharm->fill(xp, weight);
00183             _h_XpProtonCharm->fill(xp, weight);
00184             break;
00185           case PID::BQUARK:
00186             _multProton[3] += weight;
00187             _h_XpProtonBottom->fill(xp, weight);
00188             break;
00189           }
00190         }
00191       }
00192 
00193       const UnstableFinalState& ufs = applyProjection<UnstableFinalState>(e, "UFS");
00194       foreach (const Particle& p, ufs.particles()) {
00195         const double xp = p.momentum().vector3().mod()/meanBeamMom;
00196         // if in quark or antiquark hemisphere
00197         bool quark = p.momentum().vector3().dot(axis)*dot>0.;
00198         int id = abs(p.pdgId());
00199         if (id == PID::LAMBDA) {
00200           _multLambda[0] += weight;
00201           _h_XpLambdaN->fill(xp, weight);
00202           switch (flavour) {
00203           case PID::DQUARK:
00204           case PID::UQUARK:
00205           case PID::SQUARK:
00206             _multLambda[1] += weight;
00207             _h_XpLambdaLight->fill(xp, weight);
00208             if( ( quark && p.pdgId()>0 ) || ( !quark && p.pdgId()<0 ))
00209               _h_RLambda->fill(xp, weight);
00210             else
00211               _h_RLBar  ->fill(xp, weight);
00212             break;
00213           case PID::CQUARK:
00214             _multLambda[2] += weight;
00215             _h_XpLambdaCharm->fill(xp, weight);
00216             break;
00217           case PID::BQUARK:
00218             _multLambda[3] += weight;
00219             _h_XpLambdaBottom->fill(xp, weight);
00220             break;
00221           }
00222         }
00223         else if (id == 313) {
00224           _multKStar0[0] += weight;
00225           _h_XpKStar0N->fill(xp, weight);
00226           switch (flavour) {
00227           case PID::DQUARK:
00228           case PID::UQUARK:
00229           case PID::SQUARK:
00230             _multKStar0[1] += weight;
00231             _temp_XpKStar0Light->fill(xp, weight);
00232             _h_XpKStar0Light->fill(xp, weight);
00233             if ( ( quark && p.pdgId()>0 ) || ( !quark && p.pdgId()<0 ))
00234               _h_RKS0   ->fill(xp, weight);
00235             else
00236               _h_RKSBar0->fill(xp, weight);
00237             break;
00238             break;
00239           case PID::CQUARK:
00240             _multKStar0[2] += weight;
00241             _temp_XpKStar0Charm->fill(xp, weight);
00242             _h_XpKStar0Charm->fill(xp, weight);
00243             break;
00244           case PID::BQUARK:
00245             _multKStar0[3] += weight;
00246             _h_XpKStar0Bottom->fill(xp, weight);
00247             break;
00248           }
00249         }
00250         else if (id == 333) {
00251           _multPhi[0] += weight;
00252           _h_XpPhiN->fill(xp, weight);
00253           switch (flavour) {
00254           case PID::DQUARK:
00255           case PID::UQUARK:
00256           case PID::SQUARK:
00257             _multPhi[1] += weight;
00258             _h_XpPhiLight->fill(xp, weight);
00259             break;
00260           case PID::CQUARK:
00261             _multPhi[2] += weight;
00262             _h_XpPhiCharm->fill(xp, weight);
00263             break;
00264           case PID::BQUARK:
00265             _multPhi[3] += weight;
00266             _h_XpPhiBottom->fill(xp, weight);
00267             break;
00268           }
00269         }
00270         else if (id == PID::K0S || id == PID::K0L) {
00271           _multK0[0] += weight;
00272           _h_XpK0N->fill(xp, weight);
00273           switch (flavour) {
00274           case PID::DQUARK:
00275           case PID::UQUARK:
00276           case PID::SQUARK:
00277             _multK0[1] += weight;
00278             _h_XpK0Light->fill(xp, weight);
00279             break;
00280           case PID::CQUARK:
00281             _multK0[2] += weight;
00282             _h_XpK0Charm->fill(xp, weight);
00283             break;
00284           case PID::BQUARK:
00285             _multK0[3] += weight;
00286             _h_XpK0Bottom->fill(xp, weight);
00287             break;
00288           }
00289         }
00290       }
00291     }
00292 
00293 
00294     void init() {
00295       // Projections
00296       addProjection(Beam(), "Beams");
00297       addProjection(ChargedFinalState(), "FS");
00298       addProjection(UnstableFinalState(), "UFS");
00299       addProjection(InitialQuarks(), "IQF");
00300       addProjection(Thrust(FinalState()), "Thrust");
00301 
00302       _temp_XpChargedN1 = bookHisto1D( 1, 1, 1);
00303       _temp_XpChargedN2 = bookHisto1D( 2, 1, 1);
00304       _temp_XpChargedN3 = bookHisto1D( 3, 1, 1);
00305 
00306       _h_XpPiPlusN      = bookHisto1D( 1, 1, 2);
00307       _h_XpKPlusN       = bookHisto1D( 2, 1, 2);
00308       _h_XpProtonN      = bookHisto1D( 3, 1, 2);
00309       _h_XpChargedN     = bookHisto1D( 4, 1, 1);
00310       _h_XpK0N          = bookHisto1D( 5, 1, 1);
00311       _h_XpLambdaN      = bookHisto1D( 7, 1, 1);
00312       _h_XpKStar0N      = bookHisto1D( 8, 1, 1);
00313       _h_XpPhiN         = bookHisto1D( 9, 1, 1);
00314 
00315       _h_XpPiPlusLight  = bookHisto1D(10, 1, 1);
00316       _h_XpPiPlusCharm  = bookHisto1D(10, 1, 2);
00317       _h_XpPiPlusBottom = bookHisto1D(10, 1, 3);
00318       _h_XpKPlusLight   = bookHisto1D(12, 1, 1);
00319       _h_XpKPlusCharm   = bookHisto1D(12, 1, 2);
00320       _h_XpKPlusBottom  = bookHisto1D(12, 1, 3);
00321       _h_XpKStar0Light  = bookHisto1D(14, 1, 1);
00322       _h_XpKStar0Charm  = bookHisto1D(14, 1, 2);
00323       _h_XpKStar0Bottom = bookHisto1D(14, 1, 3);
00324       _h_XpProtonLight  = bookHisto1D(16, 1, 1);
00325       _h_XpProtonCharm  = bookHisto1D(16, 1, 2);
00326       _h_XpProtonBottom = bookHisto1D(16, 1, 3);
00327       _h_XpLambdaLight  = bookHisto1D(18, 1, 1);
00328       _h_XpLambdaCharm  = bookHisto1D(18, 1, 2);
00329       _h_XpLambdaBottom = bookHisto1D(18, 1, 3);
00330       _h_XpK0Light      = bookHisto1D(20, 1, 1);
00331       _h_XpK0Charm      = bookHisto1D(20, 1, 2);
00332       _h_XpK0Bottom     = bookHisto1D(20, 1, 3);
00333       _h_XpPhiLight     = bookHisto1D(22, 1, 1);
00334       _h_XpPhiCharm     = bookHisto1D(22, 1, 2);
00335       _h_XpPhiBottom    = bookHisto1D(22, 1, 3);
00336 
00337       _temp_XpKPlusCharm   = bookHisto1D("TMP/XpKPlusCharm", refData(13, 1, 1));
00338       _temp_XpKPlusLight   = bookHisto1D("TMP/XpKPlusLight", refData(13, 1, 1));
00339       _temp_XpKStar0Charm  = bookHisto1D("TMP/XpKStar0Charm", refData(15, 1, 1));
00340       _temp_XpKStar0Light  = bookHisto1D("TMP/XpKStar0Light", refData(15, 1, 1));
00341       _temp_XpProtonCharm  = bookHisto1D("TMP/XpProtonCharm", refData(17, 1, 1));
00342       _temp_XpProtonLight  = bookHisto1D("TMP/XpProtonLight", refData(17, 1, 1));
00343 
00344       _h_RPiPlus  = bookHisto1D( 26, 1, 1);
00345       _h_RPiMinus = bookHisto1D( 26, 1, 2);
00346       _h_RKS0     = bookHisto1D( 28, 1, 1);
00347       _h_RKSBar0  = bookHisto1D( 28, 1, 2);
00348       _h_RKPlus   = bookHisto1D( 30, 1, 1);
00349       _h_RKMinus  = bookHisto1D( 30, 1, 2);
00350       _h_RProton  = bookHisto1D( 32, 1, 1);
00351       _h_RPBar    = bookHisto1D( 32, 1, 2);
00352       _h_RLambda  = bookHisto1D( 34, 1, 1);
00353       _h_RLBar    = bookHisto1D( 34, 1, 2);
00354 
00355       _s_Xp_PiPl_Ch        = bookScatter2D(1, 1, 1);
00356       _s_Xp_KPl_Ch         = bookScatter2D(2, 1, 1);
00357       _s_Xp_Pr_Ch          = bookScatter2D(3, 1, 1);
00358       _s_Xp_PiPlCh_PiPlLi  = bookScatter2D(11, 1, 1);
00359       _s_Xp_PiPlBo_PiPlLi  = bookScatter2D(11, 1, 2);
00360       _s_Xp_KPlCh_KPlLi    = bookScatter2D(13, 1, 1);
00361       _s_Xp_KPlBo_KPlLi    = bookScatter2D(13, 1, 2);
00362       _s_Xp_KS0Ch_KS0Li    = bookScatter2D(15, 1, 1);
00363       _s_Xp_KS0Bo_KS0Li    = bookScatter2D(15, 1, 2);
00364       _s_Xp_PrCh_PrLi      = bookScatter2D(17, 1, 1);
00365       _s_Xp_PrBo_PrLi      = bookScatter2D(17, 1, 2);
00366       _s_Xp_LaCh_LaLi      = bookScatter2D(19, 1, 1);
00367       _s_Xp_LaBo_LaLi      = bookScatter2D(19, 1, 2);
00368       _s_Xp_K0Ch_K0Li      = bookScatter2D(21, 1, 1);
00369       _s_Xp_K0Bo_K0Li      = bookScatter2D(21, 1, 2);
00370       _s_Xp_PhiCh_PhiLi    = bookScatter2D(23, 1, 1);
00371       _s_Xp_PhiBo_PhiLi    = bookScatter2D(23, 1, 2);
00372 
00373       _s_PiM_PiP    = bookScatter2D(27, 1, 1);
00374       _s_KSBar0_KS0 = bookScatter2D(29, 1, 1);
00375       _s_KM_KP      = bookScatter2D(31, 1, 1);
00376       _s_Pr_PBar    = bookScatter2D(33, 1, 1);
00377       _s_Lam_LBar   = bookScatter2D(35, 1, 1);
00378     }
00379 
00380 
00381     /// Finalize
00382     void finalize() {
00383       // Get the ratio plots sorted out first
00384       divide(_h_XpPiPlusN, _temp_XpChargedN1, _s_Xp_PiPl_Ch);
00385       divide(_h_XpKPlusN, _temp_XpChargedN2, _s_Xp_KPl_Ch);
00386       divide(_h_XpProtonN, _temp_XpChargedN3, _s_Xp_Pr_Ch);
00387       divide(_h_XpPiPlusCharm, _h_XpPiPlusLight, _s_Xp_PiPlCh_PiPlLi);
00388       divide(_h_XpPiPlusBottom, _h_XpPiPlusLight, _s_Xp_PiPlBo_PiPlLi);
00389       divide(_temp_XpKPlusCharm , _temp_XpKPlusLight, _s_Xp_KPlCh_KPlLi);
00390       divide(_h_XpKPlusBottom, _h_XpKPlusLight, _s_Xp_KPlBo_KPlLi);
00391       divide(_temp_XpKStar0Charm, _temp_XpKStar0Light, _s_Xp_KS0Ch_KS0Li);
00392       divide(_h_XpKStar0Bottom, _h_XpKStar0Light, _s_Xp_KS0Bo_KS0Li);
00393       divide(_temp_XpProtonCharm, _temp_XpProtonLight, _s_Xp_PrCh_PrLi);
00394       divide(_h_XpProtonBottom, _h_XpProtonLight, _s_Xp_PrBo_PrLi);
00395       divide(_h_XpLambdaCharm, _h_XpLambdaLight, _s_Xp_LaCh_LaLi);
00396       divide(_h_XpLambdaBottom, _h_XpLambdaLight, _s_Xp_LaBo_LaLi);
00397       divide(_h_XpK0Charm, _h_XpK0Light, _s_Xp_K0Ch_K0Li);
00398       divide(_h_XpK0Bottom, _h_XpK0Light, _s_Xp_K0Bo_K0Li);
00399       divide(_h_XpPhiCharm, _h_XpPhiLight, _s_Xp_PhiCh_PhiLi);
00400       divide(_h_XpPhiBottom, _h_XpPhiLight, _s_Xp_PhiBo_PhiLi);
00401 
00402       // Then the leading particles
00403       divide(*_h_RPiMinus - *_h_RPiPlus, *_h_RPiMinus + *_h_RPiPlus, _s_PiM_PiP);
00404       divide(*_h_RKSBar0 - *_h_RKS0, *_h_RKSBar0 + *_h_RKS0, _s_KSBar0_KS0);
00405       divide(*_h_RKMinus - *_h_RKPlus, *_h_RKMinus + *_h_RKPlus, _s_KM_KP);
00406       divide(*_h_RProton - *_h_RPBar, *_h_RProton + *_h_RPBar, _s_Pr_PBar);
00407       divide(*_h_RLambda - *_h_RLBar, *_h_RLambda + *_h_RLBar, _s_Lam_LBar);
00408 
00409       // Then the rest
00410       scale(_h_XpPiPlusN,      1/sumOfWeights());
00411       scale(_h_XpKPlusN,       1/sumOfWeights());
00412       scale(_h_XpProtonN,      1/sumOfWeights());
00413       scale(_h_XpChargedN,     1/sumOfWeights());
00414       scale(_h_XpK0N,          1/sumOfWeights());
00415       scale(_h_XpLambdaN,      1/sumOfWeights());
00416       scale(_h_XpKStar0N,      1/sumOfWeights());
00417       scale(_h_XpPhiN,         1/sumOfWeights());
00418       scale(_h_XpPiPlusLight,  1/_SumOfudsWeights);
00419       scale(_h_XpPiPlusCharm,  1/_SumOfcWeights);
00420       scale(_h_XpPiPlusBottom, 1/_SumOfbWeights);
00421       scale(_h_XpKPlusLight,   1/_SumOfudsWeights);
00422       scale(_h_XpKPlusCharm,   1/_SumOfcWeights);
00423       scale(_h_XpKPlusBottom,  1/_SumOfbWeights);
00424       scale(_h_XpKStar0Light,  1/_SumOfudsWeights);
00425       scale(_h_XpKStar0Charm,  1/_SumOfcWeights);
00426       scale(_h_XpKStar0Bottom, 1/_SumOfbWeights);
00427       scale(_h_XpProtonLight,  1/_SumOfudsWeights);
00428       scale(_h_XpProtonCharm,  1/_SumOfcWeights);
00429       scale(_h_XpProtonBottom, 1/_SumOfbWeights);
00430       scale(_h_XpLambdaLight,  1/_SumOfudsWeights);
00431       scale(_h_XpLambdaCharm,  1/_SumOfcWeights);
00432       scale(_h_XpLambdaBottom, 1/_SumOfbWeights);
00433       scale(_h_XpK0Light,      1/_SumOfudsWeights);
00434       scale(_h_XpK0Charm,      1/_SumOfcWeights);
00435       scale(_h_XpK0Bottom,     1/_SumOfbWeights);
00436       scale(_h_XpPhiLight,     1/_SumOfudsWeights);
00437       scale(_h_XpPhiCharm ,    1/_SumOfcWeights);
00438       scale(_h_XpPhiBottom,    1/_SumOfbWeights);
00439       scale(_h_RPiPlus,        1/_SumOfudsWeights);
00440       scale(_h_RPiMinus,       1/_SumOfudsWeights);
00441       scale(_h_RKS0,           1/_SumOfudsWeights);
00442       scale(_h_RKSBar0,        1/_SumOfudsWeights);
00443       scale(_h_RKPlus,         1/_SumOfudsWeights);
00444       scale(_h_RKMinus,        1/_SumOfudsWeights);
00445       scale(_h_RProton,        1/_SumOfudsWeights);
00446       scale(_h_RPBar,          1/_SumOfudsWeights);
00447       scale(_h_RLambda,        1/_SumOfudsWeights);
00448       scale(_h_RLBar,          1/_SumOfudsWeights);
00449 
00450       // Multiplicities
00451       double avgNumPartsAll, avgNumPartsLight,avgNumPartsCharm, avgNumPartsBottom;
00452       // pi+/-
00453       // all
00454       avgNumPartsAll = _multPiPlus[0]/sumOfWeights();
00455       bookScatter2D(24, 1, 1, true)->point(0).setY(avgNumPartsAll);
00456       // light
00457       avgNumPartsLight = _multPiPlus[1]/_SumOfudsWeights;
00458       bookScatter2D(24, 1, 2, true)->point(0).setY(avgNumPartsLight);
00459       // charm
00460       avgNumPartsCharm = _multPiPlus[2]/_SumOfcWeights;
00461       bookScatter2D(24, 1, 3, true)->point(0).setY(avgNumPartsCharm);
00462       // bottom
00463       avgNumPartsBottom = _multPiPlus[3]/_SumOfbWeights;
00464       bookScatter2D(24, 1, 4, true)->point(0).setY(avgNumPartsBottom);
00465       // charm-light
00466       bookScatter2D(25, 1, 1, true)->point(0).setY(avgNumPartsCharm - avgNumPartsLight);
00467       // bottom-light
00468       bookScatter2D(25, 1, 2, true)->point(0).setY(avgNumPartsBottom - avgNumPartsLight);
00469       // K+/-
00470       // all
00471       avgNumPartsAll = _multKPlus[0]/sumOfWeights();
00472       bookScatter2D(24, 2, 1, true)->point(0).setY(avgNumPartsAll);
00473       // light
00474       avgNumPartsLight = _multKPlus[1]/_SumOfudsWeights;
00475       bookScatter2D(24, 2, 2, true)->point(0).setY(avgNumPartsLight);
00476       // charm
00477       avgNumPartsCharm = _multKPlus[2]/_SumOfcWeights;
00478       bookScatter2D(24, 2, 3, true)->point(0).setY(avgNumPartsCharm);
00479       // bottom
00480       avgNumPartsBottom = _multKPlus[3]/_SumOfbWeights;
00481       bookScatter2D(24, 2, 4, true)->point(0).setY(avgNumPartsBottom);
00482       // charm-light
00483       bookScatter2D(25, 2, 1, true)->point(0).setY(avgNumPartsCharm - avgNumPartsLight);
00484       // bottom-light
00485       bookScatter2D(25, 2, 2, true)->point(0).setY(avgNumPartsBottom - avgNumPartsLight);
00486       // K0
00487       // all
00488       avgNumPartsAll = _multK0[0]/sumOfWeights();
00489       bookScatter2D(24, 3, 1, true)->point(0).setY(avgNumPartsAll);
00490       // light
00491       avgNumPartsLight = _multK0[1]/_SumOfudsWeights;
00492       bookScatter2D(24, 3, 2, true)->point(0).setY(avgNumPartsLight);
00493       // charm
00494       avgNumPartsCharm = _multK0[2]/_SumOfcWeights;
00495       bookScatter2D(24, 3, 3, true)->point(0).setY(avgNumPartsCharm);
00496       // bottom
00497       avgNumPartsBottom = _multK0[3]/_SumOfbWeights;
00498       bookScatter2D(24, 3, 4, true)->point(0).setY(avgNumPartsBottom);
00499       // charm-light
00500       bookScatter2D(25, 3, 1, true)->point(0).setY(avgNumPartsCharm - avgNumPartsLight);
00501       // bottom-light
00502       bookScatter2D(25, 3, 2, true)->point(0).setY(avgNumPartsBottom - avgNumPartsLight);
00503       // K*0
00504       // all
00505       avgNumPartsAll = _multKStar0[0]/sumOfWeights();
00506       bookScatter2D(24, 4, 1, true)->point(0).setY(avgNumPartsAll);
00507       // light
00508       avgNumPartsLight = _multKStar0[1]/_SumOfudsWeights;
00509       bookScatter2D(24, 4, 2, true)->point(0).setY(avgNumPartsLight);
00510       // charm
00511       avgNumPartsCharm = _multKStar0[2]/_SumOfcWeights;
00512       bookScatter2D(24, 4, 3, true)->point(0).setY(avgNumPartsCharm);
00513       // bottom
00514       avgNumPartsBottom = _multKStar0[3]/_SumOfbWeights;
00515       bookScatter2D(24, 4, 4, true)->point(0).setY(avgNumPartsBottom);
00516       // charm-light
00517       bookScatter2D(25, 4, 1, true)->point(0).setY(avgNumPartsCharm - avgNumPartsLight);
00518       // bottom-light
00519       bookScatter2D(25, 4, 2, true)->point(0).setY(avgNumPartsBottom - avgNumPartsLight);
00520       // phi
00521       // all
00522       avgNumPartsAll = _multPhi[0]/sumOfWeights();
00523       bookScatter2D(24, 5, 1, true)->point(0).setY(avgNumPartsAll);
00524       // light
00525       avgNumPartsLight = _multPhi[1]/_SumOfudsWeights;
00526       bookScatter2D(24, 5, 2, true)->point(0).setY(avgNumPartsLight);
00527       // charm
00528       avgNumPartsCharm = _multPhi[2]/_SumOfcWeights;
00529       bookScatter2D(24, 5, 3, true)->point(0).setY(avgNumPartsCharm);
00530       // bottom
00531       avgNumPartsBottom = _multPhi[3]/_SumOfbWeights;
00532       bookScatter2D(24, 5, 4, true)->point(0).setY(avgNumPartsBottom);
00533       // charm-light
00534       bookScatter2D(25, 5, 1, true)->point(0).setY(avgNumPartsCharm - avgNumPartsLight);
00535       // bottom-light
00536       bookScatter2D(25, 5, 2, true)->point(0).setY(avgNumPartsBottom - avgNumPartsLight);
00537       // p
00538       // all
00539       avgNumPartsAll = _multProton[0]/sumOfWeights();
00540       bookScatter2D(24, 6, 1, true)->point(0).setY(avgNumPartsAll);
00541       // light
00542       avgNumPartsLight = _multProton[1]/_SumOfudsWeights;
00543       bookScatter2D(24, 6, 2, true)->point(0).setY(avgNumPartsLight);
00544       // charm
00545       avgNumPartsCharm = _multProton[2]/_SumOfcWeights;
00546       bookScatter2D(24, 6, 3, true)->point(0).setY(avgNumPartsCharm);
00547       // bottom
00548       avgNumPartsBottom = _multProton[3]/_SumOfbWeights;
00549       bookScatter2D(24, 6, 4, true)->point(0).setY(avgNumPartsBottom);
00550       // charm-light
00551       bookScatter2D(25, 6, 1, true)->point(0).setY(avgNumPartsCharm - avgNumPartsLight);
00552       // bottom-light
00553       bookScatter2D(25, 6, 2, true)->point(0).setY(avgNumPartsBottom - avgNumPartsLight);
00554       // Lambda
00555       // all
00556       avgNumPartsAll = _multLambda[0]/sumOfWeights();
00557       bookScatter2D(24, 7, 1, true)->point(0).setY(avgNumPartsAll);
00558       // light
00559       avgNumPartsLight = _multLambda[1]/_SumOfudsWeights;
00560       bookScatter2D(24, 7, 2, true)->point(0).setY(avgNumPartsLight);
00561       // charm
00562       avgNumPartsCharm = _multLambda[2]/_SumOfcWeights;
00563       bookScatter2D(24, 7, 3, true)->point(0).setY(avgNumPartsCharm);
00564       // bottom
00565       avgNumPartsBottom = _multLambda[3]/_SumOfbWeights;
00566       bookScatter2D(24, 7, 4, true)->point(0).setY(avgNumPartsBottom);
00567       // charm-light
00568       bookScatter2D(25, 7, 1, true)->point(0).setY(avgNumPartsCharm - avgNumPartsLight);
00569       // bottom-light
00570       bookScatter2D(25, 7, 2, true)->point(0).setY(avgNumPartsBottom - avgNumPartsLight);
00571     }
00572 
00573     //@}
00574 
00575 
00576   private:
00577 
00578     /// Store the weighted sums of numbers of charged / charged+neutral
00579     /// particles. Used to calculate average number of particles for the
00580     /// inclusive single particle distributions' normalisations.
00581     double _SumOfudsWeights, _SumOfcWeights, _SumOfbWeights;
00582     vector<double> _multPiPlus, _multKPlus, _multK0,
00583       _multKStar0, _multPhi, _multProton, _multLambda;
00584 
00585     Histo1DPtr _h_XpPiPlusSig, _h_XpPiPlusN;
00586     Histo1DPtr _h_XpKPlusSig, _h_XpKPlusN;
00587     Histo1DPtr _h_XpProtonSig, _h_XpProtonN;
00588     Histo1DPtr _h_XpChargedN;
00589     Histo1DPtr _h_XpK0N, _h_XpLambdaN;
00590     Histo1DPtr _h_XpKStar0N, _h_XpPhiN;
00591     Histo1DPtr _h_XpPiPlusLight, _h_XpPiPlusCharm, _h_XpPiPlusBottom;
00592     Histo1DPtr _h_XpKPlusLight, _h_XpKPlusCharm, _h_XpKPlusBottom;
00593     Histo1DPtr _h_XpKStar0Light, _h_XpKStar0Charm, _h_XpKStar0Bottom;
00594     Histo1DPtr _h_XpProtonLight, _h_XpProtonCharm, _h_XpProtonBottom;
00595     Histo1DPtr _h_XpLambdaLight, _h_XpLambdaCharm, _h_XpLambdaBottom;
00596     Histo1DPtr _h_XpK0Light, _h_XpK0Charm, _h_XpK0Bottom;
00597     Histo1DPtr _h_XpPhiLight, _h_XpPhiCharm, _h_XpPhiBottom;
00598 
00599     Histo1DPtr _temp_XpChargedN1, _temp_XpChargedN2, _temp_XpChargedN3;
00600     Histo1DPtr _temp_XpKPlusCharm , _temp_XpKPlusLight;
00601     Histo1DPtr _temp_XpKStar0Charm, _temp_XpKStar0Light;
00602     Histo1DPtr _temp_XpProtonCharm, _temp_XpProtonLight;
00603 
00604     Histo1DPtr _h_RPiPlus, _h_RPiMinus;
00605     Histo1DPtr _h_RKS0, _h_RKSBar0;
00606     Histo1DPtr _h_RKPlus, _h_RKMinus;
00607     Histo1DPtr _h_RProton, _h_RPBar;
00608     Histo1DPtr _h_RLambda, _h_RLBar;
00609 
00610     Scatter2DPtr _s_Xp_PiPl_Ch, _s_Xp_KPl_Ch,  _s_Xp_Pr_Ch;
00611     Scatter2DPtr _s_Xp_PiPlCh_PiPlLi, _s_Xp_PiPlBo_PiPlLi;
00612     Scatter2DPtr _s_Xp_KPlCh_KPlLi, _s_Xp_KPlBo_KPlLi;
00613     Scatter2DPtr _s_Xp_KS0Ch_KS0Li, _s_Xp_KS0Bo_KS0Li;
00614     Scatter2DPtr _s_Xp_PrCh_PrLi, _s_Xp_PrBo_PrLi;
00615     Scatter2DPtr _s_Xp_LaCh_LaLi, _s_Xp_LaBo_LaLi;
00616     Scatter2DPtr _s_Xp_K0Ch_K0Li, _s_Xp_K0Bo_K0Li;
00617     Scatter2DPtr _s_Xp_PhiCh_PhiLi, _s_Xp_PhiBo_PhiLi;
00618 
00619     Scatter2DPtr _s_PiM_PiP, _s_KSBar0_KS0, _s_KM_KP, _s_Pr_PBar, _s_Lam_LBar;
00620 
00621     //@}
00622 
00623   };
00624 
00625   // The hook for the plugin system
00626   DECLARE_RIVET_PLUGIN(SLD_1999_S3743934);
00627 
00628 }