Hemisphere.h

Go to the documentation of this file.

#ifndef PhysicsTools_Heppy_Hemisphere_h
#define PhysicsTools_Heppy_Hemisphere_h


/*  \class Hemisphere
 *
 *  Class that, given the 4-momenta of the objects in the event, 
 *  allows to split up the event into two "hemispheres" according 
 *  to different criteria (see below)/
 *
 *  Authors: Luc Pape & Filip Moortgat      Date: July 2005
 *                                          Updated: July 2006
 *                                          Updated: 15 Sept 2006
 *                                          Updated: 15 November 2006
 *                                          Updated: 01 November 2008
 *                                          Updated: 13 December 2010 (P.Nef)
 *                                          Updated: 09 February 2011
 *                                          Updated: 26 September 2011 (P.Nef)
 */


#include <vector>
#include <iostream>
#include <cmath>

namespace heppy {

class Hemisphere {

 public:

  // There are 2 constructors:
  // 1. Constructor taking as argument vectors of Px, Py, Pz and E of the objects in
  // the event that should be separated, the seeding method and the hemisphere
  // association method,
  // 2. Constructor taking as argument vectors of Px, Py, Pz and E of the objects in
  // the event that should be separated. The seeding method and the hemisphere
  // association method should then be defined by SetMethod(seeding_method, association_method).
  // 
  // Seeding method: choice of 2 inital axes
  //   1: 1st: max P ; 2nd: max P * delta R wrt first one (if delta R > 0.5)
  //   2: 2 objects who give maximal invariant mass (recommended)
  //   3: 2 objects who give maximal transverse mass
  //   4: 2 leading jets (min dR can be set by SetDRminSeed1())
  //
  // Hemisphere association method:
  //   1: maximum pt longitudinal projected on the axes
  //   2: minimal mass squared sum of the hemispheres
  //   3: minimal Lund distance (recommended)
  //
  // Note that SetMethod() also allows the seeding and/or association method to be
  // redefined for an existing hemisphere object. The GetAxis() or GetGrouping() is 
  // then recomputed using the newly defined methods.
  //
  // Some parameters possibly affecting the logic of hemisphere reconstruction can also be set
  // by the user:
  // SetNoSeed() to prevent the given object from being used as a seed
  //  (but it can be associated to hemispheres).
  // SetNoAssoc() to prevent the given object from being used in the hemisphere association
  //  (then it cannot be used as seed either).
  // ClearAllNoLists() to reset the list of NoSeed and NoAssoc objects to empty.
  // SetDRminSeed1() minimum value of DeltaR between the two hemisphere seed for seed method 1
  //   (default = 0.5). Does not affect other methods.
  // SetnItermax() maximum number of iterations allowed for the association methods (default = 100).
  //
  // Methods are provided to avoid including ISR jets into the hemispheres.
  // This is done by declaring the variable type and cut value above which an object should not
  // be included in the hemispheres.
  // Only 1 cut can be used at the time (a new declaration overwrites the previous one).
  // RejectISRPtmax() max Pt w.r.t. the hemisphere below which objects can be included
  //   (default = 10000. GeV)
  // RejectISRDRmax() max DeltaR below which objects can be included
  //   (default = 10.)

  Hemisphere(){};
  
  Hemisphere(std::vector<float> Px_vector, std::vector<float> Py_vector, std::vector<float> Pz_vector, std::vector<float> E_vector, int seed_method, int hemisphere_association_method);

  Hemisphere(std::vector<float> Px_vector, std::vector<float> Py_vector, std::vector<float> Pz_vector, std::vector<float> E_vector);

  // Destructor
  ~Hemisphere(){};

  // return Nx, Ny, Nz, P, E of the axis of group 1
  std::vector<float> getAxis1(); 
  // return Nx, Ny, Nz, P, E of the axis of group 2
  std::vector<float> getAxis2(); 

  // where Nx, Ny, Nz are the direction cosines e.g. Nx = Px/P, 
  // P is the momentum, E is the energy

  // return vector with "1" and "2"'s according to which group the object belongs 
  //   and 0 if the object is not associated to any hemisphere
  //   (order of objects in vector is same as input)
  std::vector<int> getGrouping();

  // set or overwrite the seed and association methods
  void SetMethod(int seed_method, int hemisphere_association_method) {
    seed_meth = seed_method;
    hemi_meth = hemisphere_association_method;
    status = 0;
  }

  // prevent an object from being used as a seed (but it can be associated to hemispheres)
  // (method introduced on 15/09/06)
  void SetNoSeed(int object_number) {
    Object_Noseed[object_number] = 1; 
    status = 0;
  } 

  // prevent an object from being used for hemisphere asociation (and for seeding)
  // (method introduced on 01/11/08)
  void SetNoAssoc(int object_number) {
    Object_Noassoc[object_number] = 1; 
    Object_Noseed[object_number] = 1; 
    status = 0;
  } 

  // reset the list of NoSeed and NoAssoc objects to empty
  // (method introduced on 01/11/08)
  void ClearAllNoLists() {
    for (int i = 0; i < (int)Object_Noseed.size(); ++i) {
      Object_Noassoc[i] = 0; 
      Object_Noseed[i] = 0; 
      status = 0;
    }
  }

  // set the min DeltaR for the seeds of seed method 1
  // (method introduced on 01/11/08)
  void SetDRminSeed1(float rmin) {
    dRminSeed1 = rmin;
  }

  // set the maximum number of iterations allowed for association
  // (method introduced on 01/11/08)
  void SetnItermax(int niter) {
    nItermax = niter;
  }

  // set max Pt w.r.t. the hemisphere below which objects can be included
  void RejectISRPtmax(float ptmax) {
    rejectISRPtmax = ptmax;
    rejectISRPt = 1;
    rejectISRDR = 0;
    rejectISR = 1;
    status = 0;
  }

  // set max DeltaR below which objects can be included
  void RejectISRDRmax(float drmax) {
    rejectISRDRmax = drmax;
    rejectISRDR = 1;
    rejectISRPt = 0;
    rejectISR = 1;
    status = 0;
  }

  // controls the level of debug prints
  void SetDebug(int debug) { dbg = debug; } 
  int  GetNumLoop(){return numLoop;}
 
 
 private:

  // the hemisphere separation algorithm
  int Reconstruct();
  int RejectISR(); 


  std::vector<float> Object_Px;
  std::vector<float> Object_Py;
  std::vector<float> Object_Pz;
  std::vector<float> Object_P;
  std::vector<float> Object_Pt;
  std::vector<float> Object_E;
  std::vector<float> Object_Phi;
  std::vector<float> Object_Eta;
  std::vector<int> Object_Group;
  std::vector<int> Object_Noseed;
  std::vector<int> Object_Noassoc;

  std::vector<float> Axis1;
  std::vector<float> Axis2;

  //static const float hemivsn = 1.01;
  int seed_meth;
  int hemi_meth;
  int status;
  float dRminSeed1;
  int nItermax;
  int rejectISR;
  int rejectISRPt;
  float rejectISRPtmax;
  int rejectISRDR;
  float rejectISRDRmax;
  int dbg;
  int numLoop;

};
}

#endif