Algorithm for particle flow clustering. More...

#include <PFClusterAlgo.h>

Public Types
enum	Color { NONE =0, SEED, SPECIAL }

typedef std::multimap< double, unsigned >::iterator	EH

typedef std::map< unsigned, unsigned >::const_iterator	IDH

enum	Parameter { THRESH, SEED_THRESH, PT_THRESH, SEED_PT_THRESH, CLEAN_THRESH, CLEAN_S4S1, DOUBLESPIKE_THRESH, DOUBLESPIKE_S6S2 }

typedef edm::Handle < reco::PFRecHitCollection >	PFRecHitHandle

enum	SeedState { UNKNOWN =-1, NO =0, YES =1, CLEAN =2 }

Public Member Functions
std::auto_ptr< std::vector < reco::PFCluster > > &	clusters ()

unsigned	color (unsigned rhi) const

void	doClustering (const reco::PFRecHitCollection &rechits)
	perform clustering More...

void	doClustering (const reco::PFRecHitCollection &rechits, const std::vector< bool > &mask)

void	doClustering (const PFRecHitHandle &rechitsHandle)
	perform clustering in full framework More...

void	doClustering (const PFRecHitHandle &rechitsHandle, const std::vector< bool > &mask)

void	enableDebugging (bool debug)
	set hits on which clustering will be performed More...

bool	isSeed (unsigned rhi) const

bool	masked (unsigned rhi) const

int	nNeighbours () const
	get number of neighbours for More...

double	parameter (Parameter paramtype, PFLayer::Layer layer, unsigned iCoeff=0, int iring0=0) const

	PFClusterAlgo ()
	constructor More...

int	posCalcNCrystal () const
	get number of crystals for position calculation (-1 all,5, or 9) More...

double	posCalcP1 () const
	get p1 for position calculation More...

const reco::PFRecHit &	rechit (unsigned i, const reco::PFRecHitCollection &rechits)
	More...

std::auto_ptr< std::vector < reco::PFRecHit > > &	rechitsCleaned ()

void	setCleanRBXandHPDs (bool cleanRBXandHPDs)
	Activate cleaning of HCAL RBX's and HPD's. More...

void	setHistos (TFile file, TH2F hB, TH2F *hE)
	set endcap clean threshold More...

void	setNNeighbours (int n)
	set number of neighbours for More...

void	setPosCalcNCrystal (int n)
	set number of crystals for position calculation (-1 all,5, or 9) More...

void	setPosCalcP1 (double p1)
	set p1 for position calculation More...

void	setS4S1CleanBarrel (const std::vector< double > &coeffs)

void	setS4S1CleanEndcap (const std::vector< double > &coeffs)

void	setS6S2DoubleSpikeBarrel (double cut)

void	setS6S2DoubleSpikeEndcap (double cut)

void	setShowerSigma (double sigma)
	set shower sigma for More...

void	setThreshBarrel (double thresh)
	setters -------------------------------------------------—— More...

void	setThreshCleanBarrel (double thresh)
	set barrel clean threshold More...

void	setThreshCleanEndcap (double thresh)
	set endcap clean threshold More...

void	setThreshDoubleSpikeBarrel (double thresh)
	set endcap thresholds for double spike cleaning More...

void	setThreshDoubleSpikeEndcap (double thresh)
	set endcap thresholds for double spike cleaning More...

void	setThreshEndcap (double thresh)
	set endcap threshold More...

void	setThreshPtBarrel (double thresh)

void	setThreshPtEndcap (double thresh)

void	setThreshPtSeedBarrel (double thresh)

void	setThreshPtSeedEndcap (double thresh)

void	setThreshSeedBarrel (double thresh)
	set barrel seed threshold More...

void	setThreshSeedEndcap (double thresh)
	set endcap seed threshold More...

void	setUseCornerCells (bool usecornercells)
	activate use of cells with a common corner to build topo-clusters More...

double	showerSigma () const
	get shower sigma More...

double	threshBarrel () const
	getters -------------------------------------------------—— More...

double	threshEndcap () const
	get endcap threshold More...

double	threshSeedBarrel () const
	get barrel seed threshold More...

double	threshSeedEndcap () const
	get endcap seed threshold More...

void	write ()
	write histos More...

virtual	~PFClusterAlgo ()
	destructor More...

Private Member Functions
void	buildPFClusters (const std::vector< unsigned > &cluster, const reco::PFRecHitCollection &rechits)
	build PFClusters from a topocluster More...

void	buildTopoCluster (std::vector< unsigned > &cluster, unsigned rhi, const reco::PFRecHitCollection &rechits)
	build a topocluster (recursive) More...

void	buildTopoClusters (const reco::PFRecHitCollection &rechits)
	build topoclusters around seeds More...

void	calculateClusterPosition (reco::PFCluster &cluster, reco::PFCluster &clusterwodepthcor, bool depcor=true, int posCalcNCrystal=0)
	calculate position of a cluster More...

void	cleanRBXAndHPD (const reco::PFRecHitCollection &rechits)
	Clean HCAL readout box noise and HPD discharge. More...

reco::PFRecHitRef	createRecHitRef (const reco::PFRecHitCollection &rechits, unsigned rhi)

std::pair< double, double >	dCrack (double phi, double eta)
	distance to a crack in the ECAL barrel in eta and phi direction More...

void	doClusteringWorker (const reco::PFRecHitCollection &rechits)
	perform clustering More...

void	findSeeds (const reco::PFRecHitCollection &rechits)
	look for seeds More...

void	paint (unsigned rhi, unsigned color=1)
	paint a rechit with a color. More...

Private Attributes
bool	cleanRBXandHPDs_
	option to clean HCAL RBX's and HPD's More...

std::vector< unsigned >	color_
	color, for all rechits More...

bool	debug_
	debugging on/off More...

std::multimap< double, unsigned, std::greater< double > >	eRecHits_
	indices to rechits, sorted by decreasing E (not E_T) More...

TFile *	file_

TH2F *	hBNeighbour

TH2F *	hENeighbour

std::set< unsigned >	idUsedRecHits_
	ids of rechits used in seed search More...

std::vector< bool >	mask_

std::vector< double >	minS4S1Barrel_

std::vector< double >	minS4S1Endcap_

double	minS6S2DoubleSpikeBarrel_

double	minS6S2DoubleSpikeEndcap_

int	nNeighbours_
	number of neighbours More...

std::auto_ptr< std::vector < reco::PFCluster > >	pfClusters_
	all clusters More...

std::auto_ptr< std::vector < reco::PFRecHit > >	pfRecHitsCleaned_
	particle flow rechits cleaned More...

int	posCalcNCrystal_
	number of crystals for position calculation More...

double	posCalcP1_
	parameter for position calculation More...

PFRecHitHandle	rechitsHandle_

std::vector< unsigned >	seeds_
	vector of indices for seeds. More...

std::vector< SeedState >	seedStates_
	seed state, for all rechits More...

double	showerSigma_
	sigma of shower (cm) More...

double	threshBarrel_
	barrel threshold More...

double	threshCleanBarrel_
	Barrel cleaning threshold and S4/S1 smallest fractiom. More...

double	threshCleanEndcap_
	Endcap cleaning threshold and S4/S1 smallest fractiom. More...

double	threshDoubleSpikeBarrel_
	Barrel double-spike cleaning. More...

double	threshDoubleSpikeEndcap_
	Endcap double-spike cleaning. More...

double	threshEndcap_
	endcap threshold More...

double	threshPtBarrel_

double	threshPtEndcap_

double	threshPtSeedBarrel_

double	threshPtSeedEndcap_

double	threshSeedBarrel_
	barrel seed threshold More...

double	threshSeedEndcap_
	endcap seed threshold More...

std::vector< std::vector < unsigned > >	topoClusters_
	sets of cells having one common side, and energy over threshold More...

bool	useCornerCells_
	option to use cells with a common corner to build topo-clusters More...

std::vector< bool >	usedInTopo_
	used in topo cluster? for all rechits More...

Static Private Attributes
static unsigned	prodNum_ = 1
	product number More...

Friends
std::ostream &	operator<< (std::ostream &out, const PFClusterAlgo &algo)

Detailed Description

Algorithm for particle flow clustering.

This class takes as an input a map of pointers to PFRecHit's, and creates PFCluster's from these rechits. Clustering is implemented for ECAL, HCAL, and preshower.

Author: Colin Bernet, Patrick Janot

Date: July 2006

Definition at line 31 of file PFClusterAlgo.h.

Member Typedef Documentation

typedef std::multimap<double, unsigned >::iterator PFClusterAlgo::EH

Definition at line 204 of file PFClusterAlgo.h.

typedef std::map<unsigned, unsigned >::const_iterator PFClusterAlgo::IDH

Definition at line 203 of file PFClusterAlgo.h.

typedef edm::Handle< reco::PFRecHitCollection > PFClusterAlgo::PFRecHitHandle

Definition at line 48 of file PFClusterAlgo.h.

Member Enumeration Documentation

enum PFClusterAlgo::Color

Enumerator
NONE
SEED
SPECIAL

Definition at line 154 of file PFClusterAlgo.h.

              { 
     NONE=0,
     SEED, 
     SPECIAL 
   };

enum PFClusterAlgo::Parameter

Returns: threshold, seed threshold, (gaussian width, p1 ??) for a given zone (endcap, barrel, VFCAL ??)

Enumerator
THRESH
SEED_THRESH
PT_THRESH
SEED_PT_THRESH
CLEAN_THRESH
CLEAN_S4S1
DOUBLESPIKE_THRESH
DOUBLESPIKE_S6S2

Definition at line 177 of file PFClusterAlgo.h.

                  { THRESH, 
            SEED_THRESH,
            PT_THRESH,
            SEED_PT_THRESH,
                    CLEAN_THRESH,
                    CLEAN_S4S1,
                    DOUBLESPIKE_THRESH,
                    DOUBLESPIKE_S6S2
   };

enum PFClusterAlgo::SeedState

Enumerator
UNKNOWN
NO
YES
CLEAN

Definition at line 193 of file PFClusterAlgo.h.

                  {
     UNKNOWN=-1,
     NO=0,
     YES=1,
     CLEAN=2
   };

Constructor & Destructor Documentation

PFClusterAlgo::PFClusterAlgo ( )

constructor

Definition at line 22 of file PFClusterAlgo.cc.

References file_, hBNeighbour, and hENeighbour.

                              :
   pfClusters_( new vector<reco::PFCluster> ),
   pfRecHitsCleaned_( new vector<reco::PFRecHit> ),
   threshBarrel_(0.),
   threshPtBarrel_(0.),
   threshSeedBarrel_(0.2),
   threshPtSeedBarrel_(0.),
   threshEndcap_(0.),
   threshPtEndcap_(0.),
   threshSeedEndcap_(0.6),
   threshPtSeedEndcap_(0.),
   threshCleanBarrel_(1E5),
   minS4S1Barrel_(0.),
   threshDoubleSpikeBarrel_(1E9),
   minS6S2DoubleSpikeBarrel_(-1.),
   threshCleanEndcap_(1E5),
   minS4S1Endcap_(0.),
   threshDoubleSpikeEndcap_(1E9),
   minS6S2DoubleSpikeEndcap_(-1.),
   nNeighbours_(4),
   posCalcNCrystal_(-1),
   posCalcP1_(-1),
   showerSigma_(5),
   useCornerCells_(false),
   cleanRBXandHPDs_(false),
   debug_(false) 
 {
   file_ = 0;
   hBNeighbour = 0;
   hENeighbour = 0;
 }

virtual PFClusterAlgo::~PFClusterAlgo ( )

inlinevirtual

destructor

Definition at line 39 of file PFClusterAlgo.h.

39 {;}

Member Function Documentation

void PFClusterAlgo::buildPFClusters	(	const std::vector< unsigned > &	cluster,
		const reco::PFRecHitCollection &	rechits
	)

private

build PFClusters from a topocluster

Definition at line 1002 of file PFClusterAlgo.cc.

References reco::PFCluster::addRecHitFraction(), calculateClusterPosition(), dtNoiseDBValidation_cfg::cerr, gather_cfg::cout, createRecHitRef(), debug_, delta, diffTreeTool::diff, relval_parameters_module::energy, create_public_lumi_plots::exp, cropTnPTrees::frac, i, isSeed(), reco::PFRecHit::layer(), max(), paint(), pfClusters_, posCalcNCrystal(), posCalcNCrystal_, reco::PFRecHit::position(), position, rechit(), seedStates_, showerSigma_, SPECIAL, tmp, v, and YES.

Referenced by doClusteringWorker().

 {
 
 
   //  bool debug = false;
 
 
   // several rechits may be seeds. initialize PFClusters on these seeds. 
   
   vector<reco::PFCluster> curpfclusters;
   vector<reco::PFCluster> curpfclusterswodepthcor;
   vector< unsigned > seedsintopocluster;
 
 
   for(unsigned i=0; i<topocluster.size(); i++ ) {
 
     unsigned rhi = topocluster[i];
 
     if( seedStates_[rhi] == YES ) {
 
       reco::PFCluster cluster;
       reco::PFCluster clusterwodepthcor;
 
       double fraction = 1.0; 
       
       reco::PFRecHitRef  recHitRef = createRecHitRef( rechits, rhi ); 
 
       cluster.addRecHitFraction( reco::PFRecHitFraction( recHitRef, 
                              fraction ) );
 
     // cluster.addRecHit( rhi, fraction );
 
       calculateClusterPosition( cluster,
                                 clusterwodepthcor, 
                     true );    
 
 
 //       cout<<"PFClusterAlgo: 2"<<endl;
       curpfclusters.push_back( cluster );
       curpfclusterswodepthcor.push_back( clusterwodepthcor );
 #ifdef PFLOW_DEBUG
       if(debug_) {
     cout << "PFClusterAlgo::buildPFClusters: seed "
          << rechit( rhi, rechits) <<endl;
     cout << "PFClusterAlgo::buildPFClusters: pfcluster initialized : "
          << cluster <<endl;
       }
 #endif
 
       // keep track of the seed of each topocluster
       seedsintopocluster.push_back( rhi );
     }
   }
 
   // if only one seed in the topocluster, use all crystals
   // in the position calculation (posCalcNCrystal = -1)
   // otherwise, use the user specified value
   int posCalcNCrystal = seedsintopocluster.size()>1 ? posCalcNCrystal_:-1;
   double ns2 = std::max(1.,(double)(seedsintopocluster.size())-1.);
   ns2 *= ns2;
     
   // Find iteratively the energy and position
   // of each pfcluster in the topological cluster
   unsigned iter = 0;
   unsigned niter = 50;
   double diff = ns2;
 
   // if(debug_) niter=2;
   vector<double> ener;
   vector<double> dist;
   vector<double> frac;
   vector<math::XYZVector> tmp;
 
   while ( iter++ < niter && diff > 1E-8*ns2 ) {
 
     // Store previous iteration's result and reset pfclusters     
     ener.clear();
     tmp.clear();
 
     for ( unsigned ic=0; ic<curpfclusters.size(); ic++ ) {
       ener.push_back( curpfclusters[ic].energy() );
       
       math::XYZVector v;
       v = curpfclusters[ic].position();
 
       tmp.push_back( v );
 
 #ifdef PFLOW_DEBUG
       if(debug_)  {
     cout<<"saving photon pos "<<ic<<" "<<curpfclusters[ic]<<endl;
     cout<<tmp[ic].X()<<" "<<tmp[ic].Y()<<" "<<tmp[ic].Z()<<endl;
       }
 #endif
 
       curpfclusters[ic].reset();
     }
 
     // Loop over topocluster cells
     for( unsigned irh=0; irh<topocluster.size(); irh++ ) {
       unsigned rhindex = topocluster[irh];
       
       const reco::PFRecHit& rh = rechit( rhindex, rechits);
       
       // int layer = rh.layer();
              
       dist.clear();
       frac.clear();
       double fractot = 0.;
 
       bool isaseed = isSeed(rhindex);
 
       math::XYZVector cposxyzcell;
       cposxyzcell = rh.position();
 
 #ifdef PFLOW_DEBUG
       if(debug_) { 
     cout<<rh<<endl;
     cout<<"start loop on curpfclusters"<<endl;
       }
 #endif
 
       // Loop over pfclusters
       for ( unsigned ic=0; ic<tmp.size(); ic++) {
     
 #ifdef PFLOW_DEBUG
     if(debug_) cout<<"pfcluster "<<ic<<endl;
 #endif
     
     double frc=0.;
     bool seedexclusion=true;
 
     // convert cluster coordinates to xyz
     //math::XYZVector cposxyzclust( tmp[ic].X(), tmp[ic].Y(), tmp[ic].Z() );
         // cluster position used to compute distance with cell
     math::XYZVector cposxyzclust;
         cposxyzclust = curpfclusterswodepthcor[ic].position();
 
 #ifdef PFLOW_DEBUG
     if(debug_) {
       
       cout<<"CLUSTER "<<cposxyzclust.X()<<","
           <<cposxyzclust.Y()<<","
           <<cposxyzclust.Z()<<"\t\t"
           <<"CELL "<<cposxyzcell.X()<<","
           <<cposxyzcell.Y()<<","
           <<cposxyzcell.Z()<<endl;
     }  
 #endif
     
     // Compute the distance between the current cell 
     // and the current PF cluster, normalized to a 
     // number of "sigma"
     math::XYZVector deltav = cposxyzclust;
     deltav -= cposxyzcell;
     double d = deltav.R() / showerSigma_;
     
     // if distance cell-cluster is too large, it means that 
     // we're right on the junction between 2 subdetectors (HCAL/VFCAL...)
     // in this case, distance is calculated in the xy plane
     // could also be a large supercluster... 
 #ifdef PFLOW_DEBUG
     if( d > 10. && debug_ ) { 
           paint(rhindex, SPECIAL);
       cout<<"PFClusterAlgo Warning: distance too large"<<d<<endl;
     }
 #endif
     dist.push_back( d );
 
     // the current cell is the seed from the current photon.
     if( rhindex == seedsintopocluster[ic] && seedexclusion ) {
       frc = 1.;
 #ifdef PFLOW_DEBUG
       if(debug_) cout<<"this cell is a seed for the current photon"<<endl;
 #endif
     }
     else if( isaseed && seedexclusion ) {
       frc = 0.;
 #ifdef PFLOW_DEBUG
       if(debug_) cout<<"this cell is a seed for another photon"<<endl;
 #endif
     }
     else {
       // Compute the fractions of the cell energy to be assigned to 
       // each curpfclusters in the cluster.
       frc = ener[ic] * exp ( - dist[ic]*dist[ic] / 2. );
 
 #ifdef PFLOW_DEBUG
       if(debug_) {
         cout<<"dist["<<ic<<"] "<<dist[ic]
           //        <<", sigma="<<sigma
         <<", frc="<<frc<<endl;
       }  
 #endif
     
     }
     fractot += frc;
     frac.push_back(frc);
       }      
 
       // Add the relevant fraction of the cell to the curpfclusters
 #ifdef PFLOW_DEBUG
       if(debug_) cout<<"start add cell"<<endl;
 #endif
       for ( unsigned ic=0; ic<tmp.size(); ++ic ) {
 #ifdef PFLOW_DEBUG
     if(debug_) 
       cout<<" frac["<<ic<<"] "<<frac[ic]<<" "<<fractot<<" "<<rh<<endl;
 #endif
 
     if( fractot ) 
       frac[ic] /= fractot;
     else { 
 #ifdef PFLOW_DEBUG
       if( debug_ ) {
         int layer = rh.layer();
         cerr<<"fractot = 0 ! "<<layer<<endl;
         
         for( unsigned trh=0; trh<topocluster.size(); trh++ ) {
           unsigned tindex = topocluster[trh];
           const reco::PFRecHit& rh = rechit( tindex, rechits);
           cout<<rh<<endl;
         }
 
         // assert(0)
       }
 #endif
 
       continue;
     }
 
     // if the fraction has been set to 0, the cell 
     // is now added to the cluster - careful ! (PJ, 19/07/08)
     // BUT KEEP ONLY CLOSE CELLS OTHERWISE MEMORY JUST EXPLOSES
     // (PJ, 15/09/08 <- similar to what existed before the 
         // previous bug fix, but keeps the close seeds inside, 
     // even if their fraction was set to zero.)
     // Also add a protection to keep the seed in the cluster 
     // when the latter gets far from the former. These cases
     // (about 1% of the clusters) need to be studied, as 
     // they create fake photons, in general.
     // (PJ, 16/09/08) 
         if ( dist[ic] < 10. || frac[ic] > 0.99999 ) { 
       // if ( dist[ic] > 6. ) cout << "Warning : PCluster is getting very far from its seeding cell" << endl;
       reco::PFRecHitRef  recHitRef = createRecHitRef( rechits, rhindex ); 
       reco::PFRecHitFraction rhf( recHitRef,frac[ic] );
       curpfclusters[ic].addRecHitFraction( rhf );
     }
       }
       // if(debug_) cout<<" end add cell"<<endl;
     }
 
     // Determine the new cluster position and check 
     // the distance with the previous iteration
     diff = 0.;
     for (  unsigned ic=0; ic<tmp.size(); ++ic ) {
 
       calculateClusterPosition( curpfclusters[ic], curpfclusterswodepthcor[ic], 
                                 true, posCalcNCrystal );
 #ifdef PFLOW_DEBUG
       if(debug_) cout<<"new iter "<<ic<<endl;
       if(debug_) cout<<curpfclusters[ic]<<endl;
 #endif
 
       double delta = ROOT::Math::VectorUtil::DeltaR(curpfclusters[ic].position(),tmp[ic]);
       if ( delta > diff ) diff = delta;
     }
   }
   
   // Issue a warning message if the number of iterations 
   // exceeds 50
 #ifdef PFLOW_DEBUG
   if ( iter >= 50 && debug_ ) 
     cout << "PFClusterAlgo Warning: "
      << "more than "<<niter<<" iterations in pfcluster finding: " 
      <<  setprecision(10) << diff << endl;
 #endif
   
   // There we go
   // add all clusters to the list of pfClusters.
   for(unsigned ic=0; ic<curpfclusters.size(); ic++) {
 
     calculateClusterPosition(curpfclusters[ic], curpfclusterswodepthcor[ic], 
                              true, posCalcNCrystal);
 
     pfClusters_->push_back(curpfclusters[ic]); 
   }
 }

void PFClusterAlgo::buildTopoCluster	(	std::vector< unsigned > &	cluster,
		unsigned	rhi,
		const reco::PFRecHitCollection &	rechits
	)

private

build a topocluster (recursive)

Definition at line 926 of file PFClusterAlgo.cc.

References abs, gather_cfg::cout, debug_, alignCSCRings::e, reco::PFRecHit::energy(), PFLayer::HCAL_BARREL2, i, reco::PFRecHit::layer(), masked(), reco::PFRecHit::neighbours4(), reco::PFRecHit::neighbours8(), parameter(), reco::PFRecHit::positionREP(), reco::PFRecHit::pt2(), PT_THRESH, rechit(), GOODCOLL_filter_cfg::thresh, THRESH, useCornerCells_, and usedInTopo_.

Referenced by buildTopoClusters().

                                                         {
 
 
 #ifdef PFLOW_DEBUG
   if(debug_)
     cout<<"PFClusterAlgo::buildTopoCluster in"<<endl;
 #endif
 
   const reco::PFRecHit& rh = rechit( rhi, rechits); 
 
   double e = rh.energy();
   int layer = rh.layer();
   
 
   int iring = 0;
   if (layer==PFLayer::HCAL_BARREL2 && abs(rh.positionREP().Eta())>0.34) iring= 1;
 
   double thresh = parameter( THRESH, 
                  static_cast<PFLayer::Layer>(layer), 0, iring );
   double ptThresh = parameter( PT_THRESH, 
                    static_cast<PFLayer::Layer>(layer), 0, iring );
 
 
   if( e < thresh ||  (ptThresh > 0. && rh.pt2() < ptThresh*ptThresh) ) {
 #ifdef PFLOW_DEBUG
     if(debug_)
       cout<<"return : "<<e<<"<"<<thresh<<endl; 
 #endif
     return;
   }
 
   // add hit to cluster
 
   cluster.push_back( rhi );
   // idUsedRecHits_.insert( rh.detId() );
 
   usedInTopo_[ rhi ] = true;
 
   //   cout<<" hit ptr "<<hit<<endl;
 
   // get neighbours, either with one side in common, 
   // or with one corner in common (if useCornerCells_)
   const std::vector< unsigned >& nbs 
     = useCornerCells_ ? rh.neighbours8() : rh.neighbours4();
   
   for(unsigned i=0; i<nbs.size(); i++) {
 
 //     const reco::PFRecHit& neighbour = rechit( nbs[i], rechits );
 
 //     set<unsigned>::iterator used 
 //       = idUsedRecHits_.find( neighbour.detId() );
 //     if(used != idUsedRecHits_.end() ) continue;
     
     // already used
     if( usedInTopo_[ nbs[i] ] ) {
 #ifdef PFLOW_DEBUG
       if(debug_) 
     cout<<rhi<<" used"<<endl; 
 #endif
       continue;
     }
                  
     if( !masked(nbs[i]) ) continue;
     buildTopoCluster( cluster, nbs[i], rechits );
   }
 #ifdef PFLOW_DEBUG
   if(debug_)
     cout<<"PFClusterAlgo::buildTopoCluster out"<<endl;
 #endif
 
 }

void PFClusterAlgo::buildTopoClusters ( const reco::PFRecHitCollection & rechits )

private

build topoclusters around seeds

Definition at line 882 of file PFClusterAlgo.cc.

References buildTopoCluster(), gather_cfg::cout, debug_, masked(), seeds_, topoClusters_, and usedInTopo_.

Referenced by doClusteringWorker().

                                                                             {
 
   topoClusters_.clear(); 
   
 #ifdef PFLOW_DEBUG
   if(debug_) 
     cout<<"PFClusterAlgo::buildTopoClusters start"<<endl;
 #endif
   
   for(unsigned is = 0; is<seeds_.size(); is++) {
     
     unsigned rhi = seeds_[is];
 
     if( !masked(rhi) ) continue;
     // rechit was masked to be processed
 
     // already used in a topological cluster
     if( usedInTopo_[rhi] ) {
 #ifdef PFLOW_DEBUG
       if(debug_) 
     cout<<rhi<<" used"<<endl; 
 #endif
       continue;
     }
     
     vector< unsigned > topocluster;
     buildTopoCluster( topocluster, rhi, rechits );
    
     if(topocluster.empty() ) continue;
     
     topoClusters_.push_back( topocluster );
 
   }
 
 #ifdef PFLOW_DEBUG
   if(debug_) 
     cout<<"PFClusterAlgo::buildTopoClusters done"<<endl;
 #endif
   
   return;
 }

void PFClusterAlgo::calculateClusterPosition	(	reco::PFCluster &	cluster,
		reco::PFCluster &	clusterwodepthcor,
		bool	depcor = `true`,
		int	posCalcNCrystal = `0`
	)

private

calculate position of a cluster

Definition at line 1294 of file PFClusterAlgo.cc.

References abs, dtNoiseDBValidation_cfg::cerr, gather_cfg::cout, reco::PFCluster::depthCorA_, reco::PFCluster::depthCorAp_, reco::PFCluster::depthCorB_, reco::PFCluster::depthCorBp_, reco::PFCluster::depthCorMode_, reco::PFRecHit::detId(), alignCSCRings::e, PFLayer::ECAL_BARREL, PFLayer::ECAL_ENDCAP, reco::PFRecHit::energy(), reco::CaloCluster::energy_, reco::PFRecHit::getAxisXYZ(), PFLayer::HCAL_BARREL1, PFLayer::HCAL_BARREL2, PFLayer::HCAL_ENDCAP, PFLayer::HF_EM, PFLayer::HF_HAD, reco::PFRecHit::isNeighbour4(), reco::PFRecHit::isNeighbour8(), isSeed(), reco::PFCluster::layer(), reco::PFRecHit::layer(), create_public_lumi_plots::log, max(), PFLayer::NONE, p1, posCalcNCrystal_, posCalcP1_, reco::PFRecHit::position(), reco::CaloCluster::position_, reco::PFCluster::posrep_, PFLayer::PS1, PFLayer::PS2, reco::PFCluster::rechits_, reco::PFCluster::setLayer(), mathSSE::sqrt(), threshBarrel_, threshEndcap_, x, detailsBasic3DVector::y, and detailsBasic3DVector::z.

Referenced by buildPFClusters().

                                          {
 
   if( posCalcNCrystal_ != -1 && 
       posCalcNCrystal_ != 5 && 
       posCalcNCrystal_ != 9 ) {
     throw "PFCluster::calculatePosition : posCalcNCrystal_ must be -1, 5, or 9.";
   }  
 
 
   if(!posCalcNCrystal) posCalcNCrystal = posCalcNCrystal_; 
 
   cluster.position_.SetXYZ(0,0,0);
 
   cluster.energy_ = 0;
   
   double normalize = 0;
 
   // calculate total energy, average layer, and look for seed  ---------- //
 
 
   // double layer = 0;
   map <PFLayer::Layer, double> layers; 
 
   unsigned seedIndex = 0;
   bool     seedIndexFound = false;
 
   //Colin: the following can be simplified!
 
   // loop on rechit fractions
   for (unsigned ic=0; ic<cluster.rechits_.size(); ic++ ) {
 
     unsigned rhi = cluster.rechits_[ic].recHitRef().index();
     // const reco::PFRecHit& rh = rechit( rhi, rechits );
 
     const reco::PFRecHit& rh = *(cluster.rechits_[ic].recHitRef());
     double fraction =  cluster.rechits_[ic].fraction();
 
     // Find the seed of this sub-cluster (excluding other seeds found in the topological
     // cluster, the energy fraction of which were set to 0 fpr the position determination.
     if( isSeed(rhi) && fraction > 1e-9 ) {
       seedIndex = rhi;
       seedIndexFound = true;
     }
 
     double recHitEnergy = rh.energy() * fraction;
 
     // is nan ? 
     if( recHitEnergy!=recHitEnergy ) {
       ostringstream ostr;
       edm::LogError("PFClusterAlgo")<<"rechit "<<rh.detId()<<" has a NaN energy... The input of the particle flow clustering seems to be corrupted.";
     }
 
     cluster.energy_ += recHitEnergy;
 
     // sum energy in each layer
     PFLayer::Layer layer = rh.layer();  
 
     map <PFLayer::Layer, double>:: iterator it = layers.find(layer);
 
     if (it != layers.end()) {
       it->second += recHitEnergy;
     } else {
 
       layers.insert(make_pair(layer, recHitEnergy));
 
     }
 
   }  
 
   assert(seedIndexFound);
 
   // loop over pairs to find layer with max energy          
   double Emax = 0.;
   PFLayer::Layer layer = PFLayer::NONE;
   for (map<PFLayer::Layer, double>::iterator it = layers.begin();
        it != layers.end(); ++it) {
     double e = it->second;
     if(e > Emax){ 
       Emax = e; 
       layer = it->first;
     }
   }
   
 
   //setlayer here
   cluster.setLayer( layer ); // take layer with max energy
 
   // layer /= cluster.energy_;
   // cluster.layer_ = lrintf(layer); // nearest integer
 
 
 
   double p1 =  posCalcP1_;
   if( p1 < 0 ) { 
     // automatic (and hopefully best !) determination of the parameter
     // for position determination.
     
     // Remove the ad-hoc determination of p1, and set it to the 
     // seed threshold.
     switch(cluster.layer() ) {
     case PFLayer::ECAL_BARREL:
     case PFLayer::HCAL_BARREL1:
     case PFLayer::HCAL_BARREL2:
       //    case PFLayer::HCAL_HO:
       p1 = threshBarrel_;
       break;
     case PFLayer::ECAL_ENDCAP:
     case PFLayer::HCAL_ENDCAP:
     case PFLayer::HF_EM:
     case PFLayer::HF_HAD:
     case PFLayer::PS1:
     case PFLayer::PS2:
       p1 = threshEndcap_;
       break;
 
     /*
     switch(cluster.layer() ) {
     case PFLayer::ECAL_BARREL:
       p1 = 0.004 + 0.022*cluster.energy_; // 27 feb 2006 
       break;
     case PFLayer::ECAL_ENDCAP:
       p1 = 0.014 + 0.009*cluster.energy_; // 27 feb 2006 
       break;
     case PFLayer::HCAL_BARREL1:
     case PFLayer::HCAL_BARREL2:
     case PFLayer::HCAL_ENDCAP:
     case PFLayer::HCAL_HF:
       p1 = 5.41215e-01 * log( cluster.energy_ / 1.29803e+01 );
       if(p1<0.01) p1 = 0.01;
       break;
     */
 
     default:
       cerr<<"Clusters weight_p1 -1 not yet allowed for layer "<<layer
       <<". Chose a better value in the opt file"<<endl;
       assert(0);
       break;
     }
   } 
   else if( p1< 1e-9 ) { // will divide by p1 later on
     p1 = 1e-9;
   }
   // calculate uncorrected cluster position --------------------------------
 
   reco::PFCluster::REPPoint clusterpos;   // uncorrected cluster pos 
   math::XYZPoint clusterposxyz;           // idem, xyz coord 
   math::XYZPoint firstrechitposxyz;       // pos of the rechit with highest E
 
   double maxe = -9999;
   double x = 0;
   double y = 0;
   double z = 0;
   
   for (unsigned ic=0; ic<cluster.rechits_.size(); ic++ ) {
     
     unsigned rhi = cluster.rechits_[ic].recHitRef().index();
 //     const reco::PFRecHit& rh = rechit( rhi, rechits );
 
     const reco::PFRecHit& rh = *(cluster.rechits_[ic].recHitRef());
 
     if(rhi != seedIndex) { // not the seed
       if( posCalcNCrystal == 5 ) { // pos calculated from the 5 neighbours only
     if(!rh.isNeighbour4(seedIndex) ) {
       continue;
     }
       }
       if( posCalcNCrystal == 9 ) { // pos calculated from the 9 neighbours only
     if(!rh.isNeighbour8(seedIndex) ) {
       continue;
     }
       }
     }
     double fraction =  cluster.rechits_[ic].fraction();
     double recHitEnergy = rh.energy() * fraction;
 
     double norm = fraction < 1E-9 ? 0. : max(0., log(recHitEnergy/p1 ));
     
     const math::XYZPoint& rechitposxyz = rh.position();
     
     if( recHitEnergy > maxe ) {
       firstrechitposxyz = rechitposxyz;
       maxe = recHitEnergy;
     }
 
     x += rechitposxyz.X() * norm;
     y += rechitposxyz.Y() * norm;
     z += rechitposxyz.Z() * norm;
     
     // clusterposxyz += rechitposxyz * norm;
     normalize += norm;
   }
 
   // normalize uncorrected position
   // assert(normalize);
   if( normalize < 1e-9 ) {
     //    cerr<<"--------------------"<<endl;
     //    cerr<<(*this)<<endl;
     cout << "Watch out : cluster too far from its seeding cell, set to 0,0,0" << endl;
     clusterposxyz.SetXYZ(0,0,0);
     clusterpos.SetCoordinates(0,0,0);
     return;
   }
   else {
     x /= normalize;
     y /= normalize; 
     z /= normalize; 
 
     clusterposxyz.SetCoordinates( x, y, z);
 
     clusterpos.SetCoordinates( clusterposxyz.Rho(), clusterposxyz.Eta(), clusterposxyz.Phi() );
 
   }  
 
   cluster.posrep_ = clusterpos;
 
   cluster.position_ = clusterposxyz;
 
   clusterwodepthcor = cluster;
 
 
   // correctio of the rechit position, 
   // according to the depth, only for ECAL 
 
 
   if( depcor &&   // correction requested and ECAL
       ( cluster.layer() == PFLayer::ECAL_BARREL ||       
     cluster.layer() == PFLayer::ECAL_ENDCAP ) ) {
 
 
     double corra = reco::PFCluster::depthCorA_;
     double corrb = reco::PFCluster::depthCorB_;
     if( abs(clusterpos.Eta() ) < 2.6 && 
     abs(clusterpos.Eta() ) > 1.65   ) { 
       // if crystals under preshower, correction is not the same  
       // (shower depth smaller)
       corra = reco::PFCluster::depthCorAp_;
       corrb = reco::PFCluster::depthCorBp_;
     }
 
     double depth = 0;
 
     switch( reco::PFCluster::depthCorMode_ ) {
     case 1: // for e/gamma 
       depth = corra * ( corrb + log(cluster.energy_) ); 
       break;
     case 2: // for hadrons
       depth = corra;
       break;
     default:
       cerr<<"PFClusterAlgo::calculateClusterPosition : unknown function for depth correction! "<<endl;
       assert(0);
     }
 
 
     // calculate depth vector:
     // its mag is depth
     // its direction is the cluster direction (uncorrected)
 
 //     double xdepthv = clusterposxyz.X();
 //     double ydepthv = clusterposxyz.Y();
 //     double zdepthv = clusterposxyz.Z();
 //     double mag = sqrt( xdepthv*xdepthv + 
 //             ydepthv*ydepthv + 
 //             zdepthv*zdepthv );
     
 
 //     math::XYZPoint depthv(clusterposxyz); 
 //     depthv.SetMag(depth);
     
     
     math::XYZVector depthv( clusterposxyz.X(), 
                 clusterposxyz.Y(),
                 clusterposxyz.Z() );
     depthv /= sqrt(depthv.Mag2() );
     depthv *= depth;
 
 
     // now calculate corrected cluster position:    
     math::XYZPoint clusterposxyzcor;
 
     maxe = -9999;
     x = 0;
     y = 0;
     z = 0;
     cluster.posrep_.SetXYZ(0,0,0);
     normalize = 0;
     for (unsigned ic=0; ic<cluster.rechits_.size(); ic++ ) {
 
       unsigned rhi = cluster.rechits_[ic].recHitRef().index();
 //       const reco::PFRecHit& rh = rechit( rhi, rechits );
       
       const reco::PFRecHit& rh = *(cluster.rechits_[ic].recHitRef());
 
       if(rhi != seedIndex) {
     if( posCalcNCrystal == 5 ) {
       if(!rh.isNeighbour4(seedIndex) ) {
         continue;
       }
     }
     if( posCalcNCrystal == 9 ) {
       if(!rh.isNeighbour8(seedIndex) ) {
         continue;
       }
     }
       }
 
       double fraction =  cluster.rechits_[ic].fraction();
       double recHitEnergy = rh.energy() * fraction;
       
       const math::XYZPoint&  rechitposxyz = rh.position();
 
       // rechit axis not correct ! 
       math::XYZVector rechitaxis = rh.getAxisXYZ();
       // rechitaxis -= math::XYZVector( rechitposxyz.X(), rechitposxyz.Y(), rechitposxyz.Z() );
       
       math::XYZVector rechitaxisu( rechitaxis );
       rechitaxisu /= sqrt( rechitaxis.Mag2() );
 
       math::XYZVector displacement( rechitaxisu );
       // displacement /= sqrt( displacement.Mag2() );    
       displacement *= rechitaxisu.Dot( depthv );
       
       math::XYZPoint rechitposxyzcor( rechitposxyz );
       rechitposxyzcor += displacement;
 
       if( recHitEnergy > maxe ) {
     firstrechitposxyz = rechitposxyzcor;
     maxe = recHitEnergy;
       }
 
       double norm = fraction < 1E-9 ? 0. : max(0., log(recHitEnergy/p1 ));
       
       x += rechitposxyzcor.X() * norm;
       y += rechitposxyzcor.Y() * norm;
       z += rechitposxyzcor.Z() * norm;
       
       // clusterposxyzcor += rechitposxyzcor * norm;
       normalize += norm;
     }
 
     // normalize
     if(normalize < 1e-9) {
       cerr<<"--------------------"<<endl;
       cerr<< cluster <<endl;
       assert(0);
     }
     else {
       x /= normalize;
       y /= normalize;
       z /= normalize;
       
 
       clusterposxyzcor.SetCoordinates(x,y,z);
       cluster.posrep_.SetCoordinates( clusterposxyzcor.Rho(), 
                       clusterposxyzcor.Eta(), 
                       clusterposxyzcor.Phi() );
       cluster.position_  = clusterposxyzcor;
       clusterposxyz = clusterposxyzcor;
     }
 
   }
 }

void PFClusterAlgo::cleanRBXAndHPD ( const reco::PFRecHitCollection & rechits )

private

Clean HCAL readout box noise and HPD discharge.

Definition at line 290 of file PFClusterAlgo.cc.

References abs, reco::PFRecHit::detId(), reco::PFRecHit::energy(), eRecHits_, PFLayer::HCAL_BARREL1, PFLayer::HCAL_ENDCAP, HcalDetId::ieta(), recoMuon::in, HcalDetId::iphi(), reco::PFRecHit::layer(), Association::map, mask_, masked(), reco::PFRecHit::neighbours4(), errorMatrix2Lands_multiChannel::nN, pfRecHitsCleaned_, reco::PFRecHit::position(), rechit(), mathSSE::sqrt(), and dtDQMClient_cfg::threshold.

Referenced by doClusteringWorker().

                                                                       {
 
   std::map< int, std::vector<unsigned> > hpds;
   std::map< int, std::vector<unsigned> > rbxs;
 
   for(EH ih = eRecHits_.begin(); ih != eRecHits_.end(); ih++ ) {
 
     unsigned  rhi      = ih->second; 
 
     if(! masked(rhi) ) continue;
     // rechit was asked to be processed
     const reco::PFRecHit& rhit = rechit( rhi, rechits);
     //double energy = rhit.energy();
     int layer = rhit.layer();
     if ( layer != PFLayer::HCAL_BARREL1 &&
      layer != PFLayer::HCAL_ENDCAP ) break; 
       // layer != PFLayer::HCAL_BARREL2) break; //BARREL2 for HO : need specific cleaning
     HcalDetId theHcalDetId = HcalDetId(rhit.detId());
     int ieta = theHcalDetId.ieta();
     int iphi = theHcalDetId.iphi();
     int ihpd = ieta < 0 ?  
       ( layer == PFLayer::HCAL_ENDCAP ?  -(iphi+1)/2-100 : -iphi ) : 
       ( layer == PFLayer::HCAL_ENDCAP ?   (iphi+1)/2+100 :  iphi ) ;      
     hpds[ihpd].push_back(rhi);
     int irbx = ieta < 0 ? 
       ( layer == PFLayer::HCAL_ENDCAP ?  -(iphi+5)/4 - 20 : -(iphi+5)/4 ) : 
       ( layer == PFLayer::HCAL_ENDCAP ?   (iphi+5)/4 + 20 :  (iphi+5)/4 ) ;      
     if ( irbx == 19 ) irbx = 1;
     else if ( irbx == -19 ) irbx = -1;
     else if ( irbx == 39 ) irbx = 21;
     else if ( irbx == -39 ) irbx = -21;
     rbxs[irbx].push_back(rhi);
   }
 
   // Loop on readout boxes
   for ( std::map<int, std::vector<unsigned> >::iterator itrbx = rbxs.begin();
     itrbx != rbxs.end(); ++itrbx ) { 
     if ( ( abs(itrbx->first)<20 && itrbx->second.size() > 30 ) || 
      ( abs(itrbx->first)>20 && itrbx->second.size() > 30 ) ) { 
       const std::vector<unsigned>& rhits = itrbx->second;
       double totalEta = 0.;
       double totalEtaW = 0.;
       double totalPhi = 0.;
       double totalPhiW = 0.;
       double totalEta2 = 1E-9;
       double totalEta2W = 1E-9;
       double totalPhi2 = 1E-9;
       double totalPhi2W = 1E-9;
       double totalEnergy = 0.;
       double totalEnergy2 = 1E-9;
       unsigned nSeeds = rhits.size();
       unsigned nSeeds0 = rhits.size();
       std::map< int,std::vector<unsigned> > theHPDs;
       std::multimap< double,unsigned > theEnergies;
       for ( unsigned jh=0; jh < rhits.size(); ++jh ) {
     const reco::PFRecHit& hit = rechit(rhits[jh], rechits);
     // Check if the hit is a seed
     unsigned nN = 0;
     bool isASeed = true;
     const vector<unsigned>& neighbours4 = *(& hit.neighbours4());
     for(unsigned in=0; in<neighbours4.size(); in++) {
       const reco::PFRecHit& neighbour = rechit( neighbours4[in], rechits ); 
       // one neighbour has a higher energy -> the tested rechit is not a seed
       if( neighbour.energy() > hit.energy() ) {
         --nSeeds;
         --nSeeds0;
         isASeed = false;
         break;
       } else {
         if ( neighbour.energy() > 0.4 ) ++nN;
       }
     }
     if ( isASeed && !nN ) --nSeeds0;
 
     HcalDetId theHcalDetId = HcalDetId(hit.detId());
     // int ieta = theHcalDetId.ieta();
     int iphi = theHcalDetId.iphi();
     // std::cout << "Hit : " << hit.energy() << " " << ieta << " " << iphi << std::endl;
     if ( hit.layer() == PFLayer::HCAL_BARREL1 )
       theHPDs[iphi].push_back(rhits[jh]);
     else
       theHPDs[(iphi-1)/2].push_back(rhits[jh]);
     theEnergies.insert(std::pair<double,unsigned>(hit.energy(),rhits[jh]));
     totalEnergy += hit.energy();
     totalPhi += fabs(hit.position().phi());
     totalPhiW += hit.energy()*fabs(hit.position().phi());
     totalEta += hit.position().eta();
     totalEtaW += hit.energy()*hit.position().eta();
     totalEnergy2 += hit.energy()*hit.energy();
     totalPhi2 += hit.position().phi()*hit.position().phi();
     totalPhi2W += hit.energy()*hit.position().phi()*hit.position().phi();
     totalEta2 += hit.position().eta()*hit.position().eta();
     totalEta2W += hit.energy()*hit.position().eta()*hit.position().eta();
       }
       // totalPhi /= totalEnergy;
       totalPhi /= rhits.size();
       totalEta /= rhits.size();
       totalPhiW /= totalEnergy;
       totalEtaW /= totalEnergy;
       totalPhi2 /= rhits.size();
       totalEta2 /= rhits.size();
       totalPhi2W /= totalEnergy;
       totalEta2W /= totalEnergy;
       totalPhi2 = std::sqrt(totalPhi2 - totalPhi*totalPhi);
       totalEta2 = std::sqrt(totalEta2 - totalEta*totalEta);
       totalPhi2W = std::sqrt(totalPhi2W - totalPhiW*totalPhiW);
       totalEta2W = std::sqrt(totalEta2W - totalEtaW*totalEtaW);
       totalEnergy /= rhits.size();
       totalEnergy2 /= rhits.size();
       totalEnergy2 = std::sqrt(totalEnergy2 - totalEnergy*totalEnergy);
       //if ( totalPhi2W/totalEta2W < 0.18 ) { 
       if ( nSeeds0 > 6 ) {
     unsigned nHPD15 = 0;
     for ( std::map<int, std::vector<unsigned> >::iterator itHPD = theHPDs.begin();
           itHPD != theHPDs.end(); ++itHPD ) { 
       int hpdN = itHPD->first;
       const std::vector<unsigned>& hpdHits = itHPD->second;
       if ( ( abs(hpdN) < 100 && hpdHits.size() > 14 ) || 
            ( abs(hpdN) > 100 && hpdHits.size() > 14 ) ) ++nHPD15;
     }
     if ( nHPD15 > 1 ) {
       /*
       std::cout << "Read out box numero " << itrbx->first 
             << " has " << itrbx->second.size() << " hits in it !"
             << std::endl << "sigma Eta/Phi = " << totalEta2 << " " << totalPhi2 << " " << totalPhi2/totalEta2
             << std::endl << "sigma EtaW/PhiW = " << totalEta2W << " " << totalPhi2W << " " << totalPhi2W/totalEta2W
             << std::endl << "E = " << totalEnergy << " +/- " << totalEnergy2
             << std::endl << "nSeeds = " << nSeeds << " " << nSeeds0
             << std::endl;
       for ( std::map<int, std::vector<unsigned> >::iterator itHPD = theHPDs.begin();
         itHPD != theHPDs.end(); ++itHPD ) { 
         unsigned hpdN = itHPD->first;
         const std::vector<unsigned>& hpdHits = itHPD->second;
         std::cout << "HPD number " << hpdN << " contains " << hpdHits.size() << " hits" << std::endl;
       }
       */
       std::multimap<double, unsigned >::iterator ntEn = theEnergies.end();
       --ntEn;
       unsigned nn = 0;
       double threshold = 1.;
       for ( std::multimap<double, unsigned >::iterator itEn = theEnergies.begin();
         itEn != theEnergies.end(); ++itEn ) {
         ++nn;
         if ( nn < 5 ) { 
           mask_[itEn->second] = false;
         } else if ( nn == 5 ) { 
           threshold = itEn->first * 5.;
           mask_[itEn->second] = false;
         } else { 
           if ( itEn->first < threshold ) mask_[itEn->second] = false;
         }
         if ( !masked(itEn->second) ) { 
           reco::PFRecHit theCleanedHit(rechit(itEn->second, rechits));
           //theCleanedHit.setRescale(0.);
           pfRecHitsCleaned_->push_back(theCleanedHit);
         }
         /*
         if ( !masked(itEn->second) ) 
           std::cout << "Hit Energies = " << itEn->first 
             << " " << ntEn->first/itEn->first << " masked " << std::endl;
         else 
           std::cout << "Hit Energies = " << itEn->first 
             << " " << ntEn->first/itEn->first << " kept for clustering " << std::endl;
         */
       }
     }
       }
     }
   }
 
   // Loop on hpd's
   std::map<int, std::vector<unsigned> >::iterator neighbour1;
   std::map<int, std::vector<unsigned> >::iterator neighbour2;
   std::map<int, std::vector<unsigned> >::iterator neighbour0;
   std::map<int, std::vector<unsigned> >::iterator neighbour3;
   unsigned size1 = 0;
   unsigned size2 = 0;
   for ( std::map<int, std::vector<unsigned> >::iterator ithpd = hpds.begin();
     ithpd != hpds.end(); ++ithpd ) { 
 
     const std::vector<unsigned>& rhits = ithpd->second;
     std::multimap< double,unsigned > theEnergies;
     double totalEnergy = 0.;
     double totalEnergy2 = 1E-9;
     for ( unsigned jh=0; jh < rhits.size(); ++jh ) {
       const reco::PFRecHit& hit = rechit(rhits[jh], rechits);
       totalEnergy += hit.energy();
       totalEnergy2 += hit.energy()*hit.energy();
       theEnergies.insert(std::pair<double,unsigned>(hit.energy(),rhits[jh]));
     }
     totalEnergy /= rhits.size();
     totalEnergy2 /= rhits.size();
     totalEnergy2 = std::sqrt(totalEnergy2 - totalEnergy*totalEnergy);
 
     if ( ithpd->first == 1 ) neighbour1 = hpds.find(72);
     else if ( ithpd->first == -1 ) neighbour1 = hpds.find(-72);
     else if ( ithpd->first == 101 ) neighbour1 = hpds.find(136);
     else if ( ithpd->first == -101 ) neighbour1 = hpds.find(-136);
     else neighbour1 = ithpd->first > 0 ? hpds.find(ithpd->first-1) : hpds.find(ithpd->first+1) ;
 
     if ( ithpd->first == 72 ) neighbour2 = hpds.find(1);
     else if ( ithpd->first == -72 ) neighbour2 = hpds.find(-1);
     else if ( ithpd->first == 136 ) neighbour2 = hpds.find(101);
     else if ( ithpd->first == -136 ) neighbour2 = hpds.find(-101);
     else neighbour2 = ithpd->first > 0 ? hpds.find(ithpd->first+1) : hpds.find(ithpd->first-1) ;
 
     if ( neighbour1 != hpds.end() ) { 
       if ( neighbour1->first == 1 ) neighbour0 = hpds.find(72);
       else if ( neighbour1->first == -1 ) neighbour0 = hpds.find(-72);
       else if ( neighbour1->first == 101 ) neighbour0 = hpds.find(136);
       else if ( neighbour1->first == -101 ) neighbour0 = hpds.find(-136);
       else neighbour0 = neighbour1->first > 0 ? hpds.find(neighbour1->first-1) : hpds.find(neighbour1->first+1) ;
     } 
 
     if ( neighbour2 != hpds.end() ) { 
       if ( neighbour2->first == 72 ) neighbour3 = hpds.find(1);
       else if ( neighbour2->first == -72 ) neighbour3 = hpds.find(-1);
       else if ( neighbour2->first == 136 ) neighbour3 = hpds.find(101);
       else if ( neighbour2->first == -136 ) neighbour3 = hpds.find(-101);
       else neighbour3 = neighbour2->first > 0 ? hpds.find(neighbour2->first+1) : hpds.find(neighbour2->first-1) ;
     }
 
     size1 = neighbour1 != hpds.end() ? neighbour1->second.size() : 0;
     size2 = neighbour2 != hpds.end() ? neighbour2->second.size() : 0;
 
     // Also treat the case of two neighbouring HPD's not in the same RBX
     if ( size1 > 10 ) { 
       if ( ( abs(neighbour1->first) > 100 && neighbour1->second.size() > 15 ) || 
        ( abs(neighbour1->first) < 100 && neighbour1->second.size() > 12 ) ) 
     size1 = neighbour0 != hpds.end() ? neighbour0->second.size() : 0;
     }
     if ( size2 > 10 ) { 
       if ( ( abs(neighbour2->first) > 100 && neighbour2->second.size() > 15 ) || 
        ( abs(neighbour2->first) < 100 && neighbour2->second.size() > 12 ) ) 
     size2 = neighbour3 != hpds.end() ? neighbour3->second.size() : 0;
     }
     
     if ( ( abs(ithpd->first) > 100 && ithpd->second.size() > 15 ) || 
          ( abs(ithpd->first) < 100 && ithpd->second.size() > 12 ) )
       if ( (float)(size1 + size2)/(float)ithpd->second.size() < 1.0 ) {
     /* 
     std::cout << "HPD numero " << ithpd->first 
           << " has " << ithpd->second.size() << " hits in it !" << std::endl
           << "Neighbours : " << size1 << " " << size2
           << std::endl;
     */
     std::multimap<double, unsigned >::iterator ntEn = theEnergies.end();
     --ntEn;
     unsigned nn = 0;
     double threshold = 1.;
     for ( std::multimap<double, unsigned >::iterator itEn = theEnergies.begin();
           itEn != theEnergies.end(); ++itEn ) {
       ++nn;
       if ( nn < 5 ) { 
         mask_[itEn->second] = false;
       } else if ( nn == 5 ) { 
         threshold = itEn->first * 2.5;
         mask_[itEn->second] = false;
       } else { 
         if ( itEn->first < threshold ) mask_[itEn->second] = false;
       }
       if ( !masked(itEn->second) ) { 
         reco::PFRecHit theCleanedHit(rechit(itEn->second, rechits));
         //theCleanedHit.setRescale(0.);
         pfRecHitsCleaned_->push_back(theCleanedHit);
       }
       /*
       if ( !masked(itEn->second) ) 
         std::cout << "Hit Energies = " << itEn->first 
               << " " << ntEn->first/itEn->first << " masked " << std::endl;
       else 
         std::cout << "Hit Energies = " << itEn->first 
         << " " << ntEn->first/itEn->first << " kept for clustering " << std::endl;
       */
     }
       }
   }
 
 }

std::auto_ptr< std::vector< reco::PFCluster > >& PFClusterAlgo::clusters ( )

inline

Returns: particle flow clusters

Definition at line 167 of file PFClusterAlgo.h.

References pfClusters_.

Referenced by PFRootEventManager::clustering().

168 {return pfClusters_;}

PFClusterAlgo::pfClusters_

std::auto_ptr< std::vector< reco::PFCluster > > pfClusters_

all clusters

Definition: PFClusterAlgo.h:277

unsigned PFClusterAlgo::color ( unsigned rhi ) const

Returns: color of the rechit

Definition at line 1687 of file PFClusterAlgo.cc.

References color_.

Referenced by DisplayManager::loadGRecHits(), and paint().

                                                 {
 
   if(rhi>=color_.size() ) { // rhi >= 0, since rhi is unsigned
     string err = "PFClusterAlgo::color : out of range";
     throw std::out_of_range(err);
   }
   
   return color_[rhi];
 }

reco::PFRecHitRef PFClusterAlgo::createRecHitRef	(	const reco::PFRecHitCollection &	rechits,
		unsigned	rhi
	)

private

create a reference to a rechit. in case rechitsHandle_.isValid(), this reference is permanent.

Definition at line 1722 of file PFClusterAlgo.cc.

References edm::HandleBase::isValid(), and rechitsHandle_.

Referenced by buildPFClusters().

                                {
 
   if( rechitsHandle_.isValid() ) {
     return reco::PFRecHitRef(  rechitsHandle_, rhi );
   } 
   else {
     return reco::PFRecHitRef(  &rechits, rhi );
   }
 }

std::pair< double, double > PFClusterAlgo::dCrack	(	double	phi,
		double	eta
	)

private

distance to a crack in the ECAL barrel in eta and phi direction

Definition at line 1772 of file PFClusterAlgo.cc.

References gather_cfg::cout, alignCSCRings::e, i, M_PI, min, and pi.

Referenced by findSeeds().

                                            {
 
   static double pi= M_PI;// 3.14159265358979323846;
   
   //Location of the 18 phi-cracks
   static std::vector<double> cPhi;
   if(cPhi.size()==0)
     {
       cPhi.resize(18,0);
       cPhi[0]=2.97025;
       for(unsigned i=1;i<=17;++i) cPhi[i]=cPhi[0]-2*i*pi/18;
     }
 
   //Shift of this location if eta<0
   static double delta_cPhi=0.00638;
 
   double defi; //the result
 
   //the location is shifted
   if(eta<0) phi +=delta_cPhi;
 
   if (phi>=-pi && phi<=pi){
 
     //the problem of the extrema
     if (phi<cPhi[17] || phi>=cPhi[0]){
       if (phi<0) phi+= 2*pi;
       defi = std::min(fabs(phi -cPhi[0]),fabs(phi-cPhi[17]-2*pi));          
     }
 
     //between these extrema...
     else{
       bool OK = false;
       unsigned i=16;
       while(!OK){
     if (phi<cPhi[i]){
       defi=std::min(fabs(phi-cPhi[i+1]),fabs(phi-cPhi[i]));
       OK=true;
     }
     else i-=1;
       }
     }
   }
   else{
     defi=0.;        //if there is a problem, we assum that we are in a crack
     std::cout<<"Problem in dminphi"<<std::endl;
   }
   //if(eta<0) defi=-defi;   //because of the disymetry
 
   static std::vector<double> cEta;
   if ( cEta.size() == 0 ) { 
     cEta.push_back(0.0);
     cEta.push_back(4.44747e-01)  ;   cEta.push_back(-4.44747e-01)  ;
     cEta.push_back(7.92824e-01)  ;   cEta.push_back(-7.92824e-01) ;
     cEta.push_back(1.14090e+00)  ;   cEta.push_back(-1.14090e+00)  ;
     cEta.push_back(1.47464e+00)  ;   cEta.push_back(-1.47464e+00)  ;
   }
   double deta = 999.; // the other result
 
   for ( unsigned ieta=0; ieta<cEta.size(); ++ieta ) { 
     deta = std::min(deta,fabs(eta-cEta[ieta]));
   }
   
   defi /= 0.0175;
   deta /= 0.0175;
   return std::pair<double,double>(defi,deta);
 
 }

void PFClusterAlgo::doClustering ( const reco::PFRecHitCollection & rechits )

perform clustering

Definition at line 102 of file PFClusterAlgo.cc.

References edm::HandleBase::clear(), doClusteringWorker(), mask_, rechitsHandle_, and funct::true.

Referenced by PFRootEventManager::clustering().

                                                                         {
 
   // using rechits without a Handle, clear to avoid a stale member
   rechitsHandle_.clear();
 
   // clear rechits mask
   mask_.clear();
   mask_.resize( rechits.size(), true );
 
   // perform clustering
   doClusteringWorker( rechits );
 
 }

void PFClusterAlgo::doClustering	(	const reco::PFRecHitCollection &	rechits,
		const std::vector< bool > &	mask
	)

Definition at line 116 of file PFClusterAlgo.cc.

References edm::HandleBase::clear(), doClusteringWorker(), mask_, rechitsHandle_, and funct::true.

                                                                                                       {
   // using rechits without a Handle, clear to avoid a stale member
 
   rechitsHandle_.clear();
 
   // use the specified mask, unless it doesn't match with the rechits
   mask_.clear();
 
   if (mask.size() == rechits.size()) {
       mask_.insert( mask_.end(), mask.begin(), mask.end() );
   } else {
       edm::LogError("PClusterAlgo::doClustering") << "map size should be " << rechits.size() << ". Will be reinitialized.";
       mask_.resize( rechits.size(), true );
   }
 
   // perform clustering
   doClusteringWorker( rechits );
 
 }

void PFClusterAlgo::doClustering ( const PFRecHitHandle & rechitsHandle )

perform clustering in full framework

Definition at line 66 of file PFClusterAlgo.cc.

References doClusteringWorker(), mask_, HI_PhotonSkim_cff::rechits, rechitsHandle_, and funct::true.

                                                                       {
   const reco::PFRecHitCollection& rechits = * rechitsHandle;
 
   // cache the Handle to the rechits
   rechitsHandle_ = rechitsHandle;
 
   // clear rechits mask
   mask_.clear();
   mask_.resize( rechits.size(), true );
 
   // perform clustering
   doClusteringWorker( rechits );
 }

void PFClusterAlgo::doClustering	(	const PFRecHitHandle &	rechitsHandle,
		const std::vector< bool > &	mask
	)

Definition at line 80 of file PFClusterAlgo.cc.

References doClusteringWorker(), mask_, HI_PhotonSkim_cff::rechits, rechitsHandle_, and funct::true.

                                                                                                     {
 
   const reco::PFRecHitCollection& rechits = * rechitsHandle;
 
   // cache the Handle to the rechits
   rechitsHandle_ = rechitsHandle;
 
   // use the specified mask, unless it doesn't match with the rechits
   mask_.clear();
   if (mask.size() == rechits.size()) {
       mask_.insert( mask_.end(), mask.begin(), mask.end() );
   } else {
       edm::LogError("PClusterAlgo::doClustering") << "map size should be " << rechits.size() << ". Will be reinitialized.";
       mask_.resize( rechits.size(), true );
   }
 
   // perform clustering
 
   doClusteringWorker( rechits );
 
 }

void PFClusterAlgo::doClusteringWorker ( const reco::PFRecHitCollection & rechits )

private

perform clustering

Definition at line 137 of file PFClusterAlgo.cc.

References buildPFClusters(), buildTopoClusters(), cleanRBXAndHPD(), cleanRBXandHPDs_, color_, relval_parameters_module::energy, eRecHits_, funct::false, findSeeds(), i, pfClusters_, pfRecHitsCleaned_, rechit(), seedStates_, topoClusters_, UNKNOWN, and usedInTopo_.

Referenced by doClustering().

                                                                               {
 
   if ( pfClusters_.get() )
     pfClusters_->clear();
   else 
     pfClusters_.reset( new std::vector<reco::PFCluster> );
 
   if ( pfRecHitsCleaned_.get() )
     pfRecHitsCleaned_->clear();
   else 
     pfRecHitsCleaned_.reset( new std::vector<reco::PFRecHit> );
 
   eRecHits_.clear();
 
   for ( unsigned i = 0; i < rechits.size(); i++ ){
 
     eRecHits_.insert( make_pair( rechit(i, rechits).energy(), i) );
   }
 
   color_.clear(); 
   color_.resize( rechits.size(), 0 );
 
   seedStates_.clear();
   seedStates_.resize( rechits.size(), UNKNOWN );
 
   usedInTopo_.clear();
   usedInTopo_.resize( rechits.size(), false );
 
   if ( cleanRBXandHPDs_ ) cleanRBXAndHPD( rechits);
 
   // look for seeds.
 
   findSeeds( rechits );
 
   // build topological clusters around seeds
   buildTopoClusters( rechits );
 
   // look for PFClusters inside each topological cluster (one per seed)
   
   //  int ix=0;
   //  for (reco::PFRecHitCollection::const_iterator cand =rechits.begin(); cand<rechits.end(); cand++){
   //    cout <<ix++<<" "<< cand->layer()<<endl;
   //  }
 
   for(unsigned i=0; i<topoClusters_.size(); i++) {
 
     const std::vector< unsigned >& topocluster = topoClusters_[i];
 
     buildPFClusters( topocluster, rechits ); 
 
   }
 
 }

void PFClusterAlgo::enableDebugging ( bool debug )

inline

set hits on which clustering will be performed

enable/disable debugging

Definition at line 45 of file PFClusterAlgo.h.

References debug, and debug_.

Referenced by PFRootEventManager::readOptions().

45 { debug_ = debug;}

PFClusterAlgo::debug_

bool debug_

debugging on/off

Definition: PFClusterAlgo.h:333

debug

#define debug

Definition: MEtoEDMFormat.h:34

void PFClusterAlgo::findSeeds ( const reco::PFRecHitCollection & rechits )

private

look for seeds

Definition at line 571 of file PFClusterAlgo.cc.

References abs, dtNoiseDBValidation_cfg::cerr, CLEAN, CLEAN_S4S1, CLEAN_THRESH, gather_cfg::cout, dCrack(), debug_, DOUBLESPIKE_S6S2, DOUBLESPIKE_THRESH, PFLayer::ECAL_BARREL, PFLayer::ECAL_ENDCAP, reco::PFRecHit::energy(), reco::PFRecHit::energyUp(), eRecHits_, eta(), file_, hBNeighbour, PFLayer::HCAL_BARREL1, PFLayer::HCAL_BARREL2, PFLayer::HCAL_ENDCAP, hENeighbour, PFLayer::HF_EM, PFLayer::HF_HAD, recoMuon::in, reco::PFRecHit::layer(), mask_, masked(), min, reco::PFRecHit::neighbours4(), reco::PFRecHit::neighbours8(), nNeighbours_, NO, paint(), parameter(), pfRecHitsCleaned_, phi, reco::PFRecHit::position(), reco::PFRecHit::positionREP(), PFLayer::PS1, PFLayer::PS2, reco::PFRecHit::pt2(), rechit(), SEED, SEED_PT_THRESH, SEED_THRESH, seeds_, seedStates_, and YES.

Referenced by doClusteringWorker().

                                                                      {
 
   seeds_.clear();
 
   // should replace this by the message logger.
 #ifdef PFLOW_DEBUG
   if(debug_) 
     cout<<"PFClusterAlgo::findSeeds : start"<<endl;
 #endif
 
 
   // An empty list of neighbours
   const vector<unsigned> noNeighbours(0, static_cast<unsigned>(0));
 
   // loop on rechits (sorted by decreasing energy - not E_T)
 
   for(EH ih = eRecHits_.begin(); ih != eRecHits_.end(); ih++ ) {
 
     unsigned  rhi      = ih->second; 
 
     if(! masked(rhi) ) continue;
     // rechit was asked to be processed
 
     double    rhenergy = ih->first;   
     const reco::PFRecHit& wannaBeSeed = rechit(rhi, rechits);
      
     if( seedStates_[rhi] == NO ) continue;
     // this hit was already tested, and is not a seed
  
     // determine seed energy threshold depending on the detector
     int layer = wannaBeSeed.layer();
     //for HO Ring0 and 1/2 boundary
     
     int iring = 0;
     if (layer==PFLayer::HCAL_BARREL2 && abs(wannaBeSeed.positionREP().Eta())>0.34) iring= 1;
 
     double seedThresh = parameter( SEED_THRESH, 
                    static_cast<PFLayer::Layer>(layer), 0, iring );
 
     double seedPtThresh = parameter( SEED_PT_THRESH, 
                      static_cast<PFLayer::Layer>(layer), 0., iring );
 
     double cleanThresh = parameter( CLEAN_THRESH, 
                     static_cast<PFLayer::Layer>(layer), 0, iring );
 
     double minS4S1_a = parameter( CLEAN_S4S1, 
                   static_cast<PFLayer::Layer>(layer), 0, iring );
 
     double minS4S1_b = parameter( CLEAN_S4S1, 
                   static_cast<PFLayer::Layer>(layer), 1, iring );
 
     double doubleSpikeThresh = parameter( DOUBLESPIKE_THRESH, 
                       static_cast<PFLayer::Layer>(layer), 0, iring );
 
     double doubleSpikeS6S2 = parameter( DOUBLESPIKE_S6S2, 
                     static_cast<PFLayer::Layer>(layer), 0, iring );
 
 #ifdef PFLOW_DEBUG
     if(debug_) 
       cout<<"layer:"<<layer<<" seedThresh:"<<seedThresh<<endl;
 #endif
 
 
     if( rhenergy < seedThresh || (seedPtThresh>0. && wannaBeSeed.pt2() < seedPtThresh*seedPtThresh )) {
       seedStates_[rhi] = NO; 
       continue;
     } 
 
     // Find the cell unused neighbours
     const vector<unsigned>* nbp;
     double tighterE = 1.0;
     double tighterF = 1.0;
 
     switch ( layer ) { 
     case PFLayer::ECAL_BARREL:         
     case PFLayer::ECAL_ENDCAP:       
     case PFLayer::HCAL_BARREL1:
     case PFLayer::HCAL_BARREL2:
     case PFLayer::HCAL_ENDCAP:
       tighterE = 2.0;
       tighterF = 3.0;
     case PFLayer::HF_EM:
     case PFLayer::HF_HAD:
       if( nNeighbours_ == 4 ) {
     nbp = & wannaBeSeed.neighbours4();
       }
       else if( nNeighbours_ == 8 ) {
     nbp = & wannaBeSeed.neighbours8();
       }
       else if( nNeighbours_ == 0 ) {
     nbp = & noNeighbours;
     // Allows for no clustering at all: all rechits are clusters.
     // Useful for HF
       }
       else {
     cerr<<"you're not allowed to set n neighbours to "
         <<nNeighbours_<<endl;
     assert(0);
       }
       break;
     case PFLayer::PS1:       
     case PFLayer::PS2:     
       nbp = & wannaBeSeed.neighbours4();
       break;
 
     default:
       cerr<<"CellsEF::PhotonSeeds : unknown layer "<<layer<<endl;
       assert(0);
     }
 
     const vector<unsigned>& neighbours = *nbp;
 
       
     // Select as a seed if all neighbours have a smaller energy
 
     seedStates_[rhi] = YES;
     for(unsigned in=0; in<neighbours.size(); in++) {
     
       unsigned rhj =  neighbours[in];
       // Ignore neighbours already masked
       if ( !masked(rhj) ) continue;
       const reco::PFRecHit& neighbour = rechit( rhj, rechits ); 
     
       // one neighbour has a higher energy -> the tested rechit is not a seed
       if( neighbour.energy() > wannaBeSeed.energy() ) {
     seedStates_[rhi] = NO;
     break;
       }
     }
     
     // Cleaning : check energetic, isolated seeds, likely to come from erratic noise.
     if ( file_ || wannaBeSeed.energy() > cleanThresh ) { 
       
       const vector<unsigned>& neighbours4 = *(& wannaBeSeed.neighbours4());
       // Determine the fraction of surrounding energy
       double surroundingEnergy = wannaBeSeed.energyUp();
       double neighbourEnergy = 0.;
       double layerEnergy = 0.;
       for(unsigned in4=0; in4<neighbours4.size(); in4++) {
     unsigned rhj =  neighbours4[in4];
     // Ignore neighbours already masked
     if ( !masked(rhj) ) continue;
     const reco::PFRecHit& neighbour = rechit( rhj, rechits ); 
     surroundingEnergy += neighbour.energy() + neighbour.energyUp();
     neighbourEnergy += neighbour.energy() + neighbour.energyUp();
     layerEnergy += neighbour.energy();
       }
       // Fraction 0 is the balance between EM and HAD layer for this tower
       // double fraction0 = layer == PFLayer::HF_EM || layer == PFLayer::HF_HAD ? 
       //   wannaBeSeed.energyUp()/wannaBeSeed.energy() : 1.;
       // Fraction 1 is the balance between the hit and its neighbours from both layers
       double fraction1 = surroundingEnergy/wannaBeSeed.energy();
       // Fraction 2 is the balance between the tower and the tower neighbours
       // double fraction2 = neighbourEnergy/(wannaBeSeed.energy()+wannaBeSeed.energyUp());
       // Fraction 3 is the balance between the hits and the hits neighbours in the same layer.
       // double fraction3 = layerEnergy/(wannaBeSeed.energy());
       // Mask the seed and the hit if energetic/isolated rechit
       // if ( fraction0 < minS4S1 || fraction1 < minS4S1 || fraction2 < minS4S1 || fraction3 < minS4S1 ) {
       // if ( fraction1 < minS4S1 || ( wannaBeSeed.energy() > 1.5*cleanThresh && fraction0 + fraction3 < minS4S1 ) ) {
       
       if ( file_ ) { 
     if ( layer == PFLayer::ECAL_BARREL || layer == PFLayer::HCAL_BARREL1 || layer == PFLayer::HCAL_BARREL2) { //BARREL2 for HO 
       /*
       double eta = wannaBeSeed.position().eta();
       double phi = wannaBeSeed.position().phi();
       std::pair<double,double> dcr = dCrack(phi,eta);
       double dcrmin = std::min(dcr.first, dcr.second);
       if ( dcrmin > 1. ) 
       */
       hBNeighbour->Fill(1./wannaBeSeed.energy(), fraction1);
     } else if ( fabs(wannaBeSeed.position().eta()) < 5.0  ) {
       if ( layer == PFLayer::ECAL_ENDCAP || layer == PFLayer::HCAL_ENDCAP ) { 
         if ( wannaBeSeed.energy() > 1000 ) { 
           if ( fabs(wannaBeSeed.position().eta()) < 2.8  ) 
         hENeighbour->Fill(1./wannaBeSeed.energy(), fraction1);
         }
       } else { 
         hENeighbour->Fill(1./wannaBeSeed.energy(), fraction1);
       }
     }
       }
       
       if ( wannaBeSeed.energy() > cleanThresh ) { 
     double f1Cut = minS4S1_a * log10(wannaBeSeed.energy()) + minS4S1_b;
     if ( fraction1 < f1Cut ) {
       // Double the energy cleaning threshold when close to the ECAL/HCAL - HF transition
       double eta = wannaBeSeed.position().eta();
       double phi = wannaBeSeed.position().phi();
       std::pair<double,double> dcr = dCrack(phi,eta);
       double dcrmin = layer == PFLayer::ECAL_BARREL ? std::min(dcr.first, dcr.second) : dcr.second;
       eta = fabs(eta);
       if (   eta < 5.0 &&                         // No cleaning for the HF border 
          ( ( eta < 2.85 && dcrmin > 1. ) || 
            ( rhenergy > tighterE*cleanThresh && 
              fraction1 < f1Cut/tighterF ) )  // Tighter cleaning for various cracks 
           ) { 
         seedStates_[rhi] = CLEAN;
         mask_[rhi] = false;
         reco::PFRecHit theCleanedHit(wannaBeSeed);
         //theCleanedHit.setRescale(0.);
         pfRecHitsCleaned_->push_back(theCleanedHit);
         /*
         std::cout << "A seed with E/pT/eta/phi = " << wannaBeSeed.energy() << " " << wannaBeSeed.energyUp() 
               << " " << sqrt(wannaBeSeed.pt2()) << " " << wannaBeSeed.position().eta() << " " << phi 
               << " and with surrounding fractions = " << fraction0 << " " << fraction1 
               << " " << fraction2 << " " << fraction3
               << " in layer " << layer 
               << " had been cleaned " << std::endl
               << " Distances to cracks : " << dcr.first << " " << dcr.second << std::endl
               << "(Cuts were : " << cleanThresh << " and " << f1Cut << ")" << std::endl; 
         */
       }
     }
       }
     }
 
     // Clean double spikes
     if ( mask_[rhi] && wannaBeSeed.energy() > doubleSpikeThresh ) {
       // Determine energy surrounding the seed and the most energetic neighbour
       double surroundingEnergyi = 0.;
       double enmax = -999.;
       unsigned mostEnergeticNeighbour = 0;
       const vector<unsigned>& neighbours4i = *(& wannaBeSeed.neighbours4());
       for(unsigned in4=0; in4<neighbours4i.size(); in4++) {
     unsigned rhj =  neighbours4i[in4];
     if ( !masked(rhj) ) continue;
     const reco::PFRecHit& neighbouri = rechit( rhj, rechits );
     surroundingEnergyi += neighbouri.energy();
     if ( neighbouri.energy() > enmax ) { 
       enmax = neighbouri.energy();
       mostEnergeticNeighbour = rhj;
     }
       }
       // Is there an energetic neighbour ?
       if ( enmax > 0. ) { 
     unsigned rhj = mostEnergeticNeighbour;
     const reco::PFRecHit& neighbouri = rechit( rhj, rechits );
     double surroundingEnergyj = 0.;
     //if ( mask_[rhj] && neighbouri.energy() > doubleSpikeThresh ) {
     // Determine energy surrounding the energetic neighbour
     const vector<unsigned>& neighbours4j = *(& neighbouri.neighbours4());
     for(unsigned jn4=0; jn4<neighbours4j.size(); jn4++) {
       unsigned rhk =  neighbours4j[jn4];
       const reco::PFRecHit& neighbourj = rechit( rhk, rechits ); 
       surroundingEnergyj += neighbourj.energy();
     }
     // The energy surrounding the double spike candidate 
     double surroundingEnergyFraction = 
       (surroundingEnergyi+surroundingEnergyj) / (wannaBeSeed.energy()+neighbouri.energy()) - 1.;
     if ( surroundingEnergyFraction < doubleSpikeS6S2 ) { 
       double eta = wannaBeSeed.position().eta();
       double phi = wannaBeSeed.position().phi();
       std::pair<double,double> dcr = dCrack(phi,eta);
       double dcrmin = layer == PFLayer::ECAL_BARREL ? std::min(dcr.first, dcr.second) : dcr.second;
       eta = fabs(eta);
       if (  ( eta < 5.0 && dcrmin > 1. ) ||
         ( wannaBeSeed.energy() > tighterE*doubleSpikeThresh &&
           surroundingEnergyFraction < doubleSpikeS6S2/tighterF ) ) {
         /*
         std::cout << "Double spike cleaned : Energies = " << wannaBeSeed.energy()
               << " " << neighbouri.energy() 
               << " surrounded by only " << surroundingEnergyFraction*100.
               << "% of the two spike energies at eta/phi/distance to closest crack = "  
               <<    eta << " " << phi << " " << dcrmin
               << std::endl;
         */
         // mask the seed
         seedStates_[rhi] = CLEAN;
         mask_[rhi] = false;
         reco::PFRecHit theCleanedSeed(wannaBeSeed);
         pfRecHitsCleaned_->push_back(theCleanedSeed);
         // mask the neighbour
         seedStates_[rhj] = CLEAN;
         mask_[rhj] = false;
         reco::PFRecHit theCleanedNeighbour(wannaBeSeed);
         pfRecHitsCleaned_->push_back(neighbouri);
       }
     }
       } else { 
     /*
     std::cout << "PFClusterAlgo : Double spike search : An isolated cell should have been killed " << std::endl 
           << "but is going through the double spike search!" << std::endl;
     */
       }
     }
 
     if ( seedStates_[rhi] == YES ) {
 
       // seeds_ contains the indices of all seeds. 
       seeds_.push_back( rhi );
       
       // marking the rechit
       paint(rhi, SEED);
     
       // then all neighbours cannot be seeds and are flagged as such
       for(unsigned in=0; in<neighbours.size(); in++) {
     seedStates_[ neighbours[in] ] = NO;
       }
     }
 
   }  
 
 #ifdef PFLOW_DEBUG
   if(debug_) 
     cout<<"PFClusterAlgo::findSeeds : done"<<endl;
 #endif
 }

bool PFClusterAlgo::isSeed ( unsigned rhi ) const

Returns: seed flag (not seed state) for rechit with index rhi

Definition at line 1699 of file PFClusterAlgo.cc.

References seedStates_.

Referenced by buildPFClusters(), calculateClusterPosition(), and PFRootEventManager::printRecHits().

                                              {
 
   if(rhi>=seedStates_.size() ) { // rhi >= 0, since rhi is unsigned
     string err = "PFClusterAlgo::isSeed : out of range";
     throw std::out_of_range(err);
   }
   
   return seedStates_[rhi]>0 ? true : false;
 }

bool PFClusterAlgo::masked ( unsigned rhi ) const

Returns: mask flag for the rechit

Definition at line 1676 of file PFClusterAlgo.cc.

References mask_.

Referenced by buildTopoCluster(), buildTopoClusters(), cleanRBXAndHPD(), and findSeeds().

                                              {
 
   if(rhi>=mask_.size() ) { // rhi >= 0, since rhi is unsigned
     string err = "PFClusterAlgo::masked : out of range";
     throw std::out_of_range(err);
   }
   
   return mask_[rhi];
 }

int PFClusterAlgo::nNeighbours ( ) const

inline

get number of neighbours for

Definition at line 130 of file PFClusterAlgo.h.

References nNeighbours_.

130 { return nNeighbours_;}

PFClusterAlgo::nNeighbours_

int nNeighbours_

number of neighbours

Definition: PFClusterAlgo.h:315

void PFClusterAlgo::paint	(	unsigned	rhi,
		unsigned	color = `1`
	)

private

paint a rechit with a color.

Definition at line 1710 of file PFClusterAlgo.cc.

References color(), and color_.

Referenced by buildPFClusters(), and findSeeds().

                                                        {
 
   if(rhi>=color_.size() ) { // rhi >= 0, since rhi is unsigned
     string err = "PFClusterAlgo::color : out of range";
     throw std::out_of_range(err);
   }
   
   color_[rhi] = color;
 }

double PFClusterAlgo::parameter	(	Parameter	paramtype,
		PFLayer::Layer	layer,
		unsigned	iCoeff = `0`,
		int	iring0 = `0`
	)		const

Returns: the value of a parameter of a given type, see enum Parameter, in a given layer.

Definition at line 197 of file PFClusterAlgo.cc.

References dtNoiseDBValidation_cfg::cerr, CLEAN_S4S1, CLEAN_THRESH, DOUBLESPIKE_S6S2, DOUBLESPIKE_THRESH, PFLayer::ECAL_BARREL, PFLayer::ECAL_ENDCAP, PFLayer::HCAL_BARREL1, PFLayer::HCAL_BARREL2, PFLayer::HCAL_ENDCAP, PFLayer::HF_EM, PFLayer::HF_HAD, minS4S1Barrel_, minS4S1Endcap_, minS6S2DoubleSpikeBarrel_, minS6S2DoubleSpikeEndcap_, PFLayer::PS1, PFLayer::PS2, PT_THRESH, SEED_PT_THRESH, SEED_THRESH, THRESH, threshBarrel_, threshCleanBarrel_, threshCleanEndcap_, threshDoubleSpikeBarrel_, threshDoubleSpikeEndcap_, threshEndcap_, threshPtBarrel_, threshPtEndcap_, threshPtSeedBarrel_, threshPtSeedEndcap_, threshSeedBarrel_, threshSeedEndcap_, and relativeConstraints::value.

Referenced by buildTopoCluster(), and findSeeds().

                                                      {
   
 
   double value = 0;
 
   switch( layer ) {
 
   case PFLayer::ECAL_BARREL:
   case PFLayer::HCAL_BARREL1:
   case PFLayer::HCAL_BARREL2: //  HO. 
                                  
     switch(paramtype) {
     case THRESH:
       value = (iring==0) ? threshBarrel_ : threshEndcap_; //For HO Ring0 and others
       break;
     case SEED_THRESH:
       value = (iring==0) ? threshSeedBarrel_ : threshSeedEndcap_;
       break;
     case PT_THRESH:
       value = (iring==0) ? threshPtBarrel_ : threshPtEndcap_;
       break;
     case SEED_PT_THRESH:
       value = (iring==0) ? threshPtSeedBarrel_ : threshPtSeedEndcap_;
       break;
     case CLEAN_THRESH:
       value = threshCleanBarrel_;
       break;
     case CLEAN_S4S1:
       value = minS4S1Barrel_[iCoeff];
       break;
     case DOUBLESPIKE_THRESH:
       value = threshDoubleSpikeBarrel_;
       break;
     case DOUBLESPIKE_S6S2:
       value = minS6S2DoubleSpikeBarrel_;
       break;
     default:
       cerr<<"PFClusterAlgo::parameter : unknown parameter type "
       <<paramtype<<endl;
       assert(0);
     }
     break;
   case PFLayer::ECAL_ENDCAP:
   case PFLayer::HCAL_ENDCAP:
   case PFLayer::PS1:
   case PFLayer::PS2:
   case PFLayer::HF_EM:
   case PFLayer::HF_HAD:
     // and no particle flow in VFCAL
     switch(paramtype) {
     case THRESH:
       value = threshEndcap_;
       break;
     case SEED_THRESH:
       value = threshSeedEndcap_;
       break;
     case PT_THRESH:
       value = threshPtEndcap_;
       break;
     case SEED_PT_THRESH:
       value = threshPtSeedEndcap_;
       break;
     case CLEAN_THRESH:
       value = threshCleanEndcap_;
       break;
     case CLEAN_S4S1:
       value = minS4S1Endcap_[iCoeff];
       break;
     case DOUBLESPIKE_THRESH:
       value = threshDoubleSpikeEndcap_;
       break;
     case DOUBLESPIKE_S6S2:
       value = minS6S2DoubleSpikeEndcap_;
       break;
     default:
       cerr<<"PFClusterAlgo::parameter : unknown parameter type "
       <<paramtype<<endl;
       assert(0);
     }
     break;
   default:
     cerr<<"PFClusterAlgo::parameter : unknown layer "<<layer<<endl;
     assert(0);
     break;
   }
 
   return value;
 }

int PFClusterAlgo::posCalcNCrystal ( ) const

inline

get number of crystals for position calculation (-1 all,5, or 9)

Definition at line 136 of file PFClusterAlgo.h.

References posCalcNCrystal_.

Referenced by buildPFClusters().

136 {return posCalcNCrystal_;}

PFClusterAlgo::posCalcNCrystal_

int posCalcNCrystal_

number of crystals for position calculation

Definition: PFClusterAlgo.h:318

double PFClusterAlgo::posCalcP1 ( ) const

inline

get p1 for position calculation

Definition at line 133 of file PFClusterAlgo.h.

References posCalcP1_.

133 { return posCalcP1_; }

PFClusterAlgo::posCalcP1_

double posCalcP1_

parameter for position calculation

Definition: PFClusterAlgo.h:321

const reco::PFRecHit & PFClusterAlgo::rechit	(	unsigned	i,
		const reco::PFRecHitCollection &	rechits
	)

Returns: hit with index i. performs a test on the vector range

Definition at line 1663 of file PFClusterAlgo.cc.

References i.

Referenced by buildPFClusters(), buildTopoCluster(), cleanRBXAndHPD(), doClusteringWorker(), and findSeeds().

                                                       {
 
   if(i >= rechits.size() ) { // i >= 0, since i is unsigned
     string err = "PFClusterAlgo::rechit : out of range";
     throw std::out_of_range(err);
   }
   
   return rechits[i];
 }

std::auto_ptr< std::vector< reco::PFRecHit > >& PFClusterAlgo::rechitsCleaned ( )

inline

Returns: cleaned rechits

Definition at line 171 of file PFClusterAlgo.h.

References pfRecHitsCleaned_.

172 {return pfRecHitsCleaned_;}

PFClusterAlgo::pfRecHitsCleaned_

std::auto_ptr< std::vector< reco::PFRecHit > > pfRecHitsCleaned_

particle flow rechits cleaned

Definition: PFClusterAlgo.h:280

void PFClusterAlgo::setCleanRBXandHPDs ( bool cleanRBXandHPDs )

inline

Activate cleaning of HCAL RBX's and HPD's.

Definition at line 111 of file PFClusterAlgo.h.

References cleanRBXandHPDs_.

Referenced by PFRootEventManager::readOptions().

111 { cleanRBXandHPDs_ = cleanRBXandHPDs; }

PFClusterAlgo::cleanRBXandHPDs_

bool cleanRBXandHPDs_

option to clean HCAL RBX's and HPD's

Definition: PFClusterAlgo.h:330

void PFClusterAlgo::setHistos	(	TFile *	file,
		TH2F *	hB,
		TH2F *	hE
	)

inline

set endcap clean threshold

Definition at line 93 of file PFClusterAlgo.h.

References mergeVDriftHistosByStation::file, file_, hBNeighbour, and hENeighbour.

Referenced by PFRootEventManager::readOptions().

93 {file_=file; hBNeighbour = hB; hENeighbour = hE;}

PFClusterAlgo::file_

TFile * file_

Definition: PFClusterAlgo.h:342

mergeVDriftHistosByStation.file

tuple file

Definition: mergeVDriftHistosByStation.py:54

PFClusterAlgo::hENeighbour

TH2F * hENeighbour

Definition: PFClusterAlgo.h:341

PFClusterAlgo::hBNeighbour

TH2F * hBNeighbour

Definition: PFClusterAlgo.h:340

void PFClusterAlgo::setNNeighbours ( int n )

inline

set number of neighbours for

Definition at line 96 of file PFClusterAlgo.h.

References n, and nNeighbours_.

Referenced by PFRootEventManager::readOptions().

96 { nNeighbours_ = n;}

PFClusterAlgo::nNeighbours_

int nNeighbours_

number of neighbours

Definition: PFClusterAlgo.h:315

n

int n

Definition: DTDataIntegrityTask.cc:33

void PFClusterAlgo::setPosCalcNCrystal ( int n )

inline

set number of crystals for position calculation (-1 all,5, or 9)

Definition at line 102 of file PFClusterAlgo.h.

References n, and posCalcNCrystal_.

Referenced by PFRootEventManager::readOptions().

102 { posCalcNCrystal_ = n;}

PFClusterAlgo::posCalcNCrystal_

int posCalcNCrystal_

number of crystals for position calculation

Definition: PFClusterAlgo.h:318

n

int n

Definition: DTDataIntegrityTask.cc:33

void PFClusterAlgo::setPosCalcP1 ( double p1 )

inline

set p1 for position calculation

Definition at line 99 of file PFClusterAlgo.h.

References p1, and posCalcP1_.

Referenced by PFRootEventManager::readOptions().

99 { posCalcP1_ = p1; }

p1

double p1[4]

Definition: TauolaWrapper.h:89

PFClusterAlgo::posCalcP1_

double posCalcP1_

parameter for position calculation

Definition: PFClusterAlgo.h:321

void PFClusterAlgo::setS4S1CleanBarrel ( const std::vector< double > & coeffs )

inline

Definition at line 70 of file PFClusterAlgo.h.

References minS4S1Barrel_.

Referenced by PFRootEventManager::readOptions().

70 {minS4S1Barrel_ = coeffs;}

PFClusterAlgo::minS4S1Barrel_

std::vector< double > minS4S1Barrel_

Definition: PFClusterAlgo.h:300

void PFClusterAlgo::setS4S1CleanEndcap ( const std::vector< double > & coeffs )

inline

Definition at line 86 of file PFClusterAlgo.h.

References minS4S1Endcap_.

Referenced by PFRootEventManager::readOptions().

86 {minS4S1Endcap_ = coeffs;}

PFClusterAlgo::minS4S1Endcap_

std::vector< double > minS4S1Endcap_

Definition: PFClusterAlgo.h:308

void PFClusterAlgo::setS6S2DoubleSpikeBarrel ( double cut )

inline

Definition at line 74 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::cut, and minS6S2DoubleSpikeBarrel_.

Referenced by PFRootEventManager::readOptions().

74 { minS6S2DoubleSpikeBarrel_ = cut;}

PFClusterAlgo::minS6S2DoubleSpikeBarrel_

double minS6S2DoubleSpikeBarrel_

Definition: PFClusterAlgo.h:304

GOODCOLL_filter_cfg.cut

tuple cut

Definition: GOODCOLL_filter_cfg.py:65

void PFClusterAlgo::setS6S2DoubleSpikeEndcap ( double cut )

inline

Definition at line 90 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::cut, and minS6S2DoubleSpikeEndcap_.

Referenced by PFRootEventManager::readOptions().

90 { minS6S2DoubleSpikeEndcap_ = cut;}

PFClusterAlgo::minS6S2DoubleSpikeEndcap_

double minS6S2DoubleSpikeEndcap_

Definition: PFClusterAlgo.h:312

GOODCOLL_filter_cfg.cut

tuple cut

Definition: GOODCOLL_filter_cfg.py:65

void PFClusterAlgo::setShowerSigma ( double sigma )

inline

set shower sigma for

Definition at line 105 of file PFClusterAlgo.h.

References showerSigma_.

Referenced by PFRootEventManager::readOptions().

105 { showerSigma_ = sigma;}

PFClusterAlgo::showerSigma_

double showerSigma_

sigma of shower (cm)

Definition: PFClusterAlgo.h:324

void PFClusterAlgo::setThreshBarrel ( double thresh )

inline

setters -------------------------------------------------——

set barrel threshold

Definition at line 61 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshBarrel_.

Referenced by PFRootEventManager::readOptions().

61 {threshBarrel_ = thresh;}

PFClusterAlgo::threshBarrel_

double threshBarrel_

barrel threshold

Definition: PFClusterAlgo.h:283

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshCleanBarrel ( double thresh )

inline

set barrel clean threshold

Definition at line 69 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshCleanBarrel_.

Referenced by PFRootEventManager::readOptions().

69 {threshCleanBarrel_ = thresh;}

PFClusterAlgo::threshCleanBarrel_

double threshCleanBarrel_

Barrel cleaning threshold and S4/S1 smallest fractiom.

Definition: PFClusterAlgo.h:299

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshCleanEndcap ( double thresh )

inline

set endcap clean threshold

Definition at line 85 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshCleanEndcap_.

Referenced by PFRootEventManager::readOptions().

85 {threshCleanEndcap_ = thresh;}

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

PFClusterAlgo::threshCleanEndcap_

double threshCleanEndcap_

Endcap cleaning threshold and S4/S1 smallest fractiom.

Definition: PFClusterAlgo.h:307

void PFClusterAlgo::setThreshDoubleSpikeBarrel ( double thresh )

inline

set endcap thresholds for double spike cleaning

Definition at line 73 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshDoubleSpikeBarrel_.

Referenced by PFRootEventManager::readOptions().

73 { threshDoubleSpikeBarrel_ = thresh;}

PFClusterAlgo::threshDoubleSpikeBarrel_

double threshDoubleSpikeBarrel_

Barrel double-spike cleaning.

Definition: PFClusterAlgo.h:303

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshDoubleSpikeEndcap ( double thresh )

inline

set endcap thresholds for double spike cleaning

Definition at line 89 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshDoubleSpikeEndcap_.

Referenced by PFRootEventManager::readOptions().

89 { threshDoubleSpikeEndcap_ = thresh;}

PFClusterAlgo::threshDoubleSpikeEndcap_

double threshDoubleSpikeEndcap_

Endcap double-spike cleaning.

Definition: PFClusterAlgo.h:311

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshEndcap ( double thresh )

inline

set endcap threshold

Definition at line 77 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshEndcap_.

Referenced by PFRootEventManager::readOptions().

77 {threshEndcap_ = thresh;}

PFClusterAlgo::threshEndcap_

double threshEndcap_

endcap threshold

Definition: PFClusterAlgo.h:291

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshPtBarrel ( double thresh )

inline

Definition at line 62 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshPtBarrel_.

Referenced by PFRootEventManager::readOptions().

62 {threshPtBarrel_ = thresh;}

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

PFClusterAlgo::threshPtBarrel_

double threshPtBarrel_

Definition: PFClusterAlgo.h:284

void PFClusterAlgo::setThreshPtEndcap ( double thresh )

inline

Definition at line 78 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshPtEndcap_.

Referenced by PFRootEventManager::readOptions().

78 {threshPtEndcap_ = thresh;}

PFClusterAlgo::threshPtEndcap_

double threshPtEndcap_

Definition: PFClusterAlgo.h:292

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshPtSeedBarrel ( double thresh )

inline

Definition at line 66 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshPtSeedBarrel_.

Referenced by PFRootEventManager::readOptions().

66 {threshPtSeedBarrel_ = thresh;}

PFClusterAlgo::threshPtSeedBarrel_

double threshPtSeedBarrel_

Definition: PFClusterAlgo.h:288

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshPtSeedEndcap ( double thresh )

inline

Definition at line 82 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshPtSeedEndcap_.

Referenced by PFRootEventManager::readOptions().

82 {threshPtSeedEndcap_ = thresh;}

PFClusterAlgo::threshPtSeedEndcap_

double threshPtSeedEndcap_

Definition: PFClusterAlgo.h:296

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshSeedBarrel ( double thresh )

inline

set barrel seed threshold

Definition at line 65 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshSeedBarrel_.

Referenced by PFRootEventManager::readOptions().

65 {threshSeedBarrel_ = thresh;}

PFClusterAlgo::threshSeedBarrel_

double threshSeedBarrel_

barrel seed threshold

Definition: PFClusterAlgo.h:287

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

void PFClusterAlgo::setThreshSeedEndcap ( double thresh )

inline

set endcap seed threshold

Definition at line 81 of file PFClusterAlgo.h.

References GOODCOLL_filter_cfg::thresh, and threshSeedEndcap_.

Referenced by PFRootEventManager::readOptions().

81 {threshSeedEndcap_ = thresh;}

GOODCOLL_filter_cfg.thresh

tuple thresh

Definition: GOODCOLL_filter_cfg.py:74

PFClusterAlgo::threshSeedEndcap_

double threshSeedEndcap_

endcap seed threshold

Definition: PFClusterAlgo.h:295

void PFClusterAlgo::setUseCornerCells ( bool usecornercells )

inline

activate use of cells with a common corner to build topo-clusters

Definition at line 108 of file PFClusterAlgo.h.

References useCornerCells_.

Referenced by PFRootEventManager::readOptions().

108 { useCornerCells_ = usecornercells;}

PFClusterAlgo::useCornerCells_

bool useCornerCells_

option to use cells with a common corner to build topo-clusters

Definition: PFClusterAlgo.h:327

double PFClusterAlgo::showerSigma ( ) const

inline

get shower sigma

Definition at line 139 of file PFClusterAlgo.h.

References showerSigma_.

139 { return showerSigma_ ;}

PFClusterAlgo::showerSigma_

double showerSigma_

sigma of shower (cm)

Definition: PFClusterAlgo.h:324

double PFClusterAlgo::threshBarrel ( ) const

inline

getters -------------------------------------------------——

get barrel threshold

Definition at line 116 of file PFClusterAlgo.h.

References threshBarrel_.

Referenced by DisplayManager::createGRecHit().

116 {return threshBarrel_;}

PFClusterAlgo::threshBarrel_

double threshBarrel_

barrel threshold

Definition: PFClusterAlgo.h:283

double PFClusterAlgo::threshEndcap ( ) const

inline

get endcap threshold

Definition at line 123 of file PFClusterAlgo.h.

References threshEndcap_.

Referenced by DisplayManager::createGRecHit().

123 {return threshEndcap_;}

PFClusterAlgo::threshEndcap_

double threshEndcap_

endcap threshold

Definition: PFClusterAlgo.h:291

double PFClusterAlgo::threshSeedBarrel ( ) const

inline

get barrel seed threshold

Definition at line 119 of file PFClusterAlgo.h.

References threshSeedBarrel_.

119 {return threshSeedBarrel_;}

PFClusterAlgo::threshSeedBarrel_

double threshSeedBarrel_

barrel seed threshold

Definition: PFClusterAlgo.h:287

double PFClusterAlgo::threshSeedEndcap ( ) const

inline

get endcap seed threshold

Definition at line 126 of file PFClusterAlgo.h.

References threshSeedEndcap_.

126 {return threshSeedEndcap_;}

PFClusterAlgo::threshSeedEndcap_

double threshSeedEndcap_

endcap seed threshold

Definition: PFClusterAlgo.h:295

void PFClusterAlgo::write ( void )

write histos

Definition at line 55 of file PFClusterAlgo.cc.

References gather_cfg::cout, and file_.

Referenced by PFRootEventManager::write().

                      { 
 
   if ( file_ ) { 
     file_->Write();
     cout << "Benchmark output written to file " << file_->GetName() << endl;
     file_->Close();
   }
 
 }

Friends And Related Function Documentation

std::ostream& operator<<	(	std::ostream &	out,
		const PFClusterAlgo &	algo
	)

friend

Member Data Documentation

bool PFClusterAlgo::cleanRBXandHPDs_

private

option to clean HCAL RBX's and HPD's

Definition at line 330 of file PFClusterAlgo.h.

Referenced by doClusteringWorker(), and setCleanRBXandHPDs().

std::vector< unsigned > PFClusterAlgo::color_

private

color, for all rechits

Definition at line 259 of file PFClusterAlgo.h.

Referenced by color(), doClusteringWorker(), and paint().

bool PFClusterAlgo::debug_

private

debugging on/off

Definition at line 333 of file PFClusterAlgo.h.

Referenced by buildPFClusters(), buildTopoCluster(), buildTopoClusters(), enableDebugging(), and findSeeds().

std::multimap<double, unsigned, std::greater<double> > PFClusterAlgo::eRecHits_

private

indices to rechits, sorted by decreasing E (not E_T)

Definition at line 252 of file PFClusterAlgo.h.

Referenced by cleanRBXAndHPD(), doClusteringWorker(), and findSeeds().

TFile* PFClusterAlgo::file_

private

Definition at line 342 of file PFClusterAlgo.h.

Referenced by findSeeds(), PFClusterAlgo(), setHistos(), and write().

TH2F* PFClusterAlgo::hBNeighbour

private

Definition at line 340 of file PFClusterAlgo.h.

Referenced by findSeeds(), PFClusterAlgo(), and setHistos().

TH2F* PFClusterAlgo::hENeighbour

private

Definition at line 341 of file PFClusterAlgo.h.

Referenced by findSeeds(), PFClusterAlgo(), and setHistos().

std::set<unsigned> PFClusterAlgo::idUsedRecHits_

private

ids of rechits used in seed search

Definition at line 249 of file PFClusterAlgo.h.

std::vector< bool > PFClusterAlgo::mask_

private

mask, for all rechits. only masked rechits will be clustered by default, all rechits are masked. see setMask function.

Definition at line 256 of file PFClusterAlgo.h.

Referenced by cleanRBXAndHPD(), doClustering(), findSeeds(), and masked().

std::vector<double> PFClusterAlgo::minS4S1Barrel_

private

Definition at line 300 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), and setS4S1CleanBarrel().

std::vector<double> PFClusterAlgo::minS4S1Endcap_

private

Definition at line 308 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), and setS4S1CleanEndcap().

double PFClusterAlgo::minS6S2DoubleSpikeBarrel_

private

Definition at line 304 of file PFClusterAlgo.h.

Referenced by parameter(), and setS6S2DoubleSpikeBarrel().

double PFClusterAlgo::minS6S2DoubleSpikeEndcap_

private

Definition at line 312 of file PFClusterAlgo.h.

Referenced by parameter(), and setS6S2DoubleSpikeEndcap().

int PFClusterAlgo::nNeighbours_

private

number of neighbours

Definition at line 315 of file PFClusterAlgo.h.

Referenced by findSeeds(), nNeighbours(), operator<<(), and setNNeighbours().

std::auto_ptr< std::vector<reco::PFCluster> > PFClusterAlgo::pfClusters_

private

all clusters

particle flow clusters

Definition at line 277 of file PFClusterAlgo.h.

Referenced by buildPFClusters(), clusters(), doClusteringWorker(), and operator<<().

std::auto_ptr< std::vector<reco::PFRecHit> > PFClusterAlgo::pfRecHitsCleaned_

private

particle flow rechits cleaned

Definition at line 280 of file PFClusterAlgo.h.

Referenced by cleanRBXAndHPD(), doClusteringWorker(), findSeeds(), and rechitsCleaned().

int PFClusterAlgo::posCalcNCrystal_

private

number of crystals for position calculation

Definition at line 318 of file PFClusterAlgo.h.

Referenced by buildPFClusters(), calculateClusterPosition(), operator<<(), posCalcNCrystal(), and setPosCalcNCrystal().

double PFClusterAlgo::posCalcP1_

private

parameter for position calculation

Definition at line 321 of file PFClusterAlgo.h.

Referenced by calculateClusterPosition(), operator<<(), posCalcP1(), and setPosCalcP1().

unsigned PFClusterAlgo::prodNum_ = 1

staticprivate

product number

Definition at line 337 of file PFClusterAlgo.h.

PFRecHitHandle PFClusterAlgo::rechitsHandle_

private

Definition at line 246 of file PFClusterAlgo.h.

Referenced by createRecHitRef(), and doClustering().

std::vector< unsigned > PFClusterAlgo::seeds_

private

vector of indices for seeds.

Definition at line 268 of file PFClusterAlgo.h.

Referenced by buildTopoClusters(), and findSeeds().

std::vector< SeedState > PFClusterAlgo::seedStates_

private

seed state, for all rechits

Definition at line 262 of file PFClusterAlgo.h.

Referenced by buildPFClusters(), doClusteringWorker(), findSeeds(), and isSeed().

double PFClusterAlgo::showerSigma_

private

sigma of shower (cm)

Definition at line 324 of file PFClusterAlgo.h.

Referenced by buildPFClusters(), operator<<(), setShowerSigma(), and showerSigma().

double PFClusterAlgo::threshBarrel_

private

barrel threshold

Definition at line 283 of file PFClusterAlgo.h.

Referenced by calculateClusterPosition(), operator<<(), parameter(), setThreshBarrel(), and threshBarrel().

double PFClusterAlgo::threshCleanBarrel_

private

Barrel cleaning threshold and S4/S1 smallest fractiom.

Definition at line 299 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), and setThreshCleanBarrel().

double PFClusterAlgo::threshCleanEndcap_

private

Endcap cleaning threshold and S4/S1 smallest fractiom.

Definition at line 307 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), and setThreshCleanEndcap().

double PFClusterAlgo::threshDoubleSpikeBarrel_

private

Barrel double-spike cleaning.

Definition at line 303 of file PFClusterAlgo.h.

Referenced by parameter(), and setThreshDoubleSpikeBarrel().

double PFClusterAlgo::threshDoubleSpikeEndcap_

private

Endcap double-spike cleaning.

Definition at line 311 of file PFClusterAlgo.h.

Referenced by parameter(), and setThreshDoubleSpikeEndcap().

double PFClusterAlgo::threshEndcap_

private

endcap threshold

Definition at line 291 of file PFClusterAlgo.h.

Referenced by calculateClusterPosition(), operator<<(), parameter(), setThreshEndcap(), and threshEndcap().

double PFClusterAlgo::threshPtBarrel_

private

Definition at line 284 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), and setThreshPtBarrel().

double PFClusterAlgo::threshPtEndcap_

private

Definition at line 292 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), and setThreshPtEndcap().

double PFClusterAlgo::threshPtSeedBarrel_

private

Definition at line 288 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), and setThreshPtSeedBarrel().

double PFClusterAlgo::threshPtSeedEndcap_

private

Definition at line 296 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), and setThreshPtSeedEndcap().

double PFClusterAlgo::threshSeedBarrel_

private

barrel seed threshold

Definition at line 287 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), setThreshSeedBarrel(), and threshSeedBarrel().

double PFClusterAlgo::threshSeedEndcap_

private

endcap seed threshold

Definition at line 295 of file PFClusterAlgo.h.

Referenced by operator<<(), parameter(), setThreshSeedEndcap(), and threshSeedEndcap().

std::vector< std::vector< unsigned > > PFClusterAlgo::topoClusters_

private

sets of cells having one common side, and energy over threshold

Definition at line 271 of file PFClusterAlgo.h.

Referenced by buildTopoClusters(), and doClusteringWorker().

bool PFClusterAlgo::useCornerCells_

private

option to use cells with a common corner to build topo-clusters

Definition at line 327 of file PFClusterAlgo.h.

Referenced by buildTopoCluster(), operator<<(), and setUseCornerCells().

std::vector< bool > PFClusterAlgo::usedInTopo_

private

used in topo cluster? for all rechits

Definition at line 265 of file PFClusterAlgo.h.

Referenced by buildTopoCluster(), buildTopoClusters(), and doClusteringWorker().

Public Types

Public Member Functions

Private Member Functions

Private Attributes

Static Private Attributes

Friends

Detailed Description

Member Typedef Documentation

Member Enumeration Documentation

Constructor & Destructor Documentation

Member Function Documentation

Friends And Related Function Documentation

Member Data Documentation