baseAnalysis.h

Go to the documentation of this file.

#ifndef baseAnalysis_H
#define baseAnalysis_H

#include "HepMC/GenEvent.h"

#include <iostream>
#include <stdio.h>
#include "HepMC/IO_GenEvent.h"
#include "HepMC/GenEvent.h"
#include "HepMC/IO_AsciiParticles.h"
#include "HepMC/SimpleVector.h"

#include "HepPDT/ParticleDataTable.hh"

#include "fastjet/PseudoJet.hh"
#include "fastjet/ClusterSequence.hh"
#include "fastjet/JetDefinition.hh"
#include "fastjet/SISConePlugin.hh"

#include "TH1.h"

using namespace std;

enum RETURN{FAILURE=-1, SUCCESS=1, CHARGED=2, NOTCHARGED=3};

/**
@class baseAnalysis.h
@brief Base class for all  HepMCanalysis classes to implement common 
       functionality and interfaces.

@author Cano Ay Dec 2008 */


class baseAnalysis
{
 public:
  baseAnalysis();
  virtual ~baseAnalysis();
  
  int Finalize();
  int Finalize(TFile* output);
  
  virtual int Init(double maxeta=2.5, double minpt=0.5) { std::cout << "baseAnalysis: WARNING: You have called the dummy function Init()." << std::endl; return 0;};
  virtual int Process(HepMC::GenEvent* hepmcevt) { std::cout << "baseAnalysis: WARNING: You have called the dummy function Process()." << std::endl; return 0;} ;

  /** Initialize FastJet Algorithm*/
  int InitJetFinder(double coneRadius, double overlapThreshold, double jet_ptmin);
  int FindJetableParticles(HepMC::GenEvent* hepmcevt);
  int FindJet(HepMC::GenEvent* hepmcevt);
  int DeleteJetObject(HepMC::GenEvent* hepmcevt);

  /** Set the Output filename*/  
  inline int setOutpuFileName(const char* filename){ m_outputFileName=filename; return 0;};
  
  /** Set the directory name in Output root file*/  
  inline int setOutpuRootDir(const char* dirname){ m_outputRootDir=dirname; return 0;};

  /** Check some special neutral Particles */
  inline bool IsNeutrino(int pid) { if( abs(pid)==12 || abs(pid)==14 || abs(pid)==16) return true; else return false; };
  inline bool IsGamma(int pid) { if( abs(pid)==22 ) return true; else return false; };
  inline bool IsNeutron(int pid) { if( abs(pid)==2112 ) return true; else return false; };
  inline bool IsK0(int pid) { if( abs(pid)==130 || abs(pid)==310 || abs(pid)==311 ) return true; else return false; };
  inline bool IsPi0(int pid) { if( abs(pid)==111 || abs(pid)==113 || abs(pid)==221 ) return true; else return false; };
    
  /** Check if neutral particle*/  
  inline int chargedParticle(int pid) { 
      if( IsNeutrino(pid) || IsGamma(pid) || IsNeutron(pid) || IsK0(pid) || IsPi0(pid) || abs(pid)== 94 )
      {return NOTCHARGED;} 
    else 
      {return CHARGED;}
    return 0;};
  
  /** Set the maximum allowed eta range*/  
  inline int setMaxEtaCut(double maxeta){ m_max_eta=maxeta; return 0;};
  
  /** Set maximum pt of tracks */  
  inline int setMinPtCut(double minpt){ m_min_pt=minpt; return 0;};
  
  /** Initialization of histograms  */ 
  TH1D * initHist(string name, string title, string xlabel, int nrBins=100, double xmin=0., double xmax=100.) ;
  TH1D * initHistVariableBin(string name, string title, string xlabel, int nbin, Double_t nbinRange[]) ;
  
  /** check if mother decayed into daugther */ 
  bool checkDaughter(HepMC::GenParticle* mother, HepMC::GenParticle* daughter, int maxGenerations=-1) ;
  
  /** check if the track comes from a specific particle e.g. pid(W-Boson)=24 */ 
  bool trackfromPID(int pid,HepMC::GenParticle* track, int maxGenerations=-1) ;

  /** calculate the rapidity of a particle */
  inline double getRapidity(HepMC::GenEvent::particle_const_iterator p) {
    double e = (*p)->momentum().e();
    double pz = (*p)->momentum().pz();
    double rapidity = 0.5 * log((e + pz) / (e - pz));
    return rapidity;
  };

  
 protected:
  
  /** jet finding */
  fastjet::JetDefinition::Plugin*   m_plugin;
  fastjet::JetDefinition*           m_jetDef;

  std::vector<fastjet::PseudoJet>  m_input_particles;
  std::vector<fastjet::PseudoJet> m_inclusive_jets;
  fastjet::ClusterSequence* m_clust_seq;

  vector< TH1D * > m_histVector;
  vector< TH1D * > m_histVectorVariableBin;  
  string m_outputFileName;
  string m_outputRootDir;

  //for jet finder
  double m_coneRadius;
  double m_overlapThreshold;
  double m_jet_ptmin;
  bool m_Jetfinder_enabled;
  
  // specify the maximum eta of tracks
  double m_max_eta;
  double m_min_pt;
  
  const HepPDT::ParticleDataTable* m_particleTable;
};

#endif

---