rivet is hosted by Hepforge, IPPP Durham
TASSO_1990_S2148048.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/Thrust.hh"
00005 #include "Rivet/Projections/Sphericity.hh"
00006 #include "Rivet/Projections/ChargedFinalState.hh"
00007 /// @todo Include more projections as required, e.g. ChargedFinalState, FastJets, ZFinder...
00008 
00009 namespace Rivet {
00010 
00011 
00012   class TASSO_1990_S2148048 : public Analysis {
00013   public:
00014 
00015     /// @name Constructors etc.
00016     //@{
00017 
00018     /// Constructor
00019     TASSO_1990_S2148048()
00020       : Analysis("TASSO_1990_S2148048")
00021     {
00022       /// @todo Set whether your finalize method needs the generator cross section
00023       //_sumWPassed = 0;
00024     }
00025 
00026     //@}
00027 
00028 
00029   public:
00030 
00031     /// @name Analysis methods
00032     //@{
00033 
00034     /// Book histograms and initialise projections before the run
00035     void init() {
00036       const ChargedFinalState cfs(-MAXRAPIDITY, MAXRAPIDITY, 0.1/GeV);
00037       addProjection(cfs, "CFS");
00038 
00039       //// Beams -- needed for x_p calculation
00040       //addProjection(Beam(), "Beams");
00041 
00042       // Thrust
00043       addProjection(Thrust(cfs), "Thrust");
00044 
00045       // For Sphericity and the like
00046       addProjection(Sphericity(cfs), "Sphericity");
00047 
00048       // Histos
00049       int offset = 0;
00050       switch (int(sqrtS()/GeV)) {
00051         case 14:
00052           offset = 0;
00053           break;
00054         case 22:
00055           offset = 1;
00056           break;
00057         case 35:
00058           offset = 2;
00059           break;
00060         case 44:
00061           offset = 3;
00062           break;
00063       }
00064       //_h_xp         = bookHisto1D( 2, 1, 1+offset);
00065       _h_sphericity = bookHisto1D( 6, 1, 1+offset);
00066       _h_aplanarity = bookHisto1D( 7, 1, 1+offset);
00067       _h_thrust     = bookHisto1D( 8, 1, 1+offset);
00068     }
00069 
00070 
00071 
00072     /// Perform the per-event analysis
00073     void analyze(const Event& event) {
00074       const double weight = event.weight();
00075       const ChargedFinalState& cfs = applyProjection<ChargedFinalState>(event, "CFS");
00076 
00077       //// Get beams and average beam momentum
00078       //const ParticlePair& beams = applyProjection<Beam>(event, "Beams").beams();
00079       //const double meanBeamMom = ( beams.first.momentum().vector3().mod() +
00080                                    //beams.second.momentum().vector3().mod() ) / 2.0;
00081 
00082       // TASSO hadronic event selection TODO: move this into a trigger definition
00083       // See page 2 in publication
00084       // Condition 1)  --- require at least 5 (4) 'good' tracks
00085       int nch = cfs.particles().size();
00086       if ( (int(sqrtS()/GeV) > 27 && nch < 5) || (int(sqrtS()/GeV) <= 27 && nch < 4 ) ) {
00087         MSG_DEBUG("Failed # good tracks cut: " << nch);
00088         vetoEvent;
00089       }
00090       // Condition 2) ---
00091       // Condition 5) --- scalar momentum (not pT!!!) sum >= 0.265*s
00092       double momsum = 0.0;
00093       foreach (const Particle& p, cfs.particles()) {
00094         const double mom = p.momentum().vector3().mod();
00095         momsum += mom;
00096       }
00097       if (momsum <=0.265 * sqrtS()/GeV) {
00098         MSG_DEBUG("Failed pTsum cut: " << momsum << " < " << 0.265 * sqrtS()/GeV);
00099         vetoEvent;
00100       }
00101 
00102       // Raise counter for events that pass trigger conditions
00103       //_sumWPassed += event.weight();
00104 
00105       const Thrust& thrust = applyProjection<Thrust>(event, "Thrust");
00106       //const Vector3 & thrustAxis = thrust.thrustAxis ();
00107       //double theta = thrustAxis.theta();
00108       //if ( fabs(cos(theta)) >= 0.8 ) {
00109         //MSG_DEBUG("Failed thrust angle cut: " << fabs(cos(theta)));
00110         //vetoEvent;
00111       //}
00112 
00113       const Sphericity& sphericity = applyProjection<Sphericity>(event, "Sphericity");
00114 
00115       //// Fill histograms in order of appearance in paper
00116       //foreach (const Particle& p, cfs.particles()) {
00117         //// Get momentum and energy of each particle.
00118         //const Vector3 mom3 = p.momentum().vector3();
00119         //// Scaled momenta.
00120         //const double mom = mom3.mod();
00121         //const double scaledMom = mom/meanBeamMom;
00122         //_h_xp->fill(scaledMom, weight);
00123       //}
00124       //
00125       _h_sphericity->fill(sphericity.sphericity(), weight);
00126       _h_aplanarity->fill(sphericity.aplanarity(), weight);
00127       _h_thrust->fill(thrust.thrust(), weight);
00128     }
00129 
00130 
00131     /// Normalise histograms etc., after the run
00132     void finalize() {
00133       //scale(_h_xp, _sumWPassed/(crossSection()*sumOfWeights()));
00134       normalize(_h_sphericity);
00135       normalize(_h_aplanarity);
00136       normalize(_h_thrust    );
00137     }
00138 
00139     //@}
00140 
00141 
00142   private:
00143 
00144     // Data members like post-cuts event weight counters go here
00145 
00146     //double _sumWPassed;
00147 
00148   private:
00149 
00150     /// @name Histograms
00151     //@{
00152     Histo1DPtr _h_xp        ;
00153     Histo1DPtr _h_sphericity;
00154     Histo1DPtr _h_aplanarity;
00155     Histo1DPtr _h_thrust    ;
00156     //@}
00157 
00158 
00159   };
00160 
00161 
00162 
00163   // The hook for the plugin system
00164   DECLARE_RIVET_PLUGIN(TASSO_1990_S2148048);
00165 
00166 }