TreeFiller/src/FillerSuperClusters.cc

// $Id: FillerSuperClusters.cc,v 1.13 2010/09/19 23:47:53 bendavid Exp $

#include "MitProd/TreeFiller/interface/FillerSuperClusters.h"
#include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
#include "DataFormats/EgammaReco/interface/SuperCluster.h"
#include "MitAna/DataTree/interface/BasicCluster.h"
#include "MitAna/DataTree/interface/SuperClusterCol.h"
#include "MitAna/DataTree/interface/Names.h"
#include "MitProd/ObjectService/interface/ObjectService.h"
#include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaTowerIsolation.h"
#include "HiggsAnalysis/HiggsToGammaGamma/interface/PhotonFix.h"

using namespace std;
using namespace edm;
using namespace mithep;

//--------------------------------------------------------------------------------------------------
FillerSuperClusters::FillerSuperClusters(const ParameterSet &cfg, const char *name, bool active) : 
  BaseFiller(cfg,name,active),
  edmName_(Conf().getUntrackedParameter<string>("edmName","hybridSuperClusters")),
  mitName_(Conf().getUntrackedParameter<string>("mitName","SuperClusters")),
  basicClusterMapName_(Conf().getUntrackedParameter<string>("basicClusterMapName", 
                                                            "BasicClusterMap")),
  caloTowerDetIdMapName_(Conf().getUntrackedParameter<string>("caloTowerDetIdMapName", 
                                                            "CaloTowerDetIdMap")),                                                            
  superClusterMapName_(Conf().getUntrackedParameter<string>("superClusterMapName", 
                                                            "SuperClusterMap")),
  superClusterIdMapName_(Conf().getUntrackedParameter<string>("superClusterIdMapName", 
                                                              "SuperClusterIdMap")),
  caloTowerName_(Conf().getUntrackedParameter<string>("caloTowerName","towerMaker")),
  basicClusterMap_(0),
  caloTowerDetIdMap_(0),
  superClusters_(new mithep::SuperClusterArr(25)),
  superClusterMap_(new mithep::SuperClusterMap),
  superClusterIdMap_(new mithep::SuperClusterIdMap)
{
  // Constructor.
}

//--------------------------------------------------------------------------------------------------
FillerSuperClusters::~FillerSuperClusters()
{
  // Destructor.

  delete superClusters_;
  delete superClusterMap_;
  delete superClusterIdMap_;
}

//--------------------------------------------------------------------------------------------------
void FillerSuperClusters::BookDataBlock(TreeWriter &tws)
{
  // Add super cluster branch to tree and get pointers to maps.

  tws.AddBranch(mitName_,&superClusters_);
  OS()->add<SuperClusterArr>(superClusters_,mitName_);

  if (!basicClusterMapName_.empty()) {
    basicClusterMap_ = OS()->get<BasicClusterMap>(basicClusterMapName_);
    if (basicClusterMap_)
      AddBranchDep(mitName_,basicClusterMap_->GetBrName());
  }
  
  if (!caloTowerDetIdMapName_.empty()) {
    caloTowerDetIdMap_ = OS()->get<CaloTowerDetIdMap>(caloTowerDetIdMapName_);
    if (caloTowerDetIdMap_)
      AddBranchDep(mitName_,caloTowerDetIdMap_->GetBrName());
  }
  
  if (!superClusterMapName_.empty()) {
    superClusterMap_->SetBrName(mitName_);
    OS()->add<SuperClusterMap>(superClusterMap_,superClusterMapName_);
  }
  if (!superClusterIdMapName_.empty()) {
    superClusterIdMap_->SetBrName(mitName_);
    OS()->add<SuperClusterIdMap>(superClusterIdMap_,superClusterIdMapName_);
  }
}

//--------------------------------------------------------------------------------------------------
void FillerSuperClusters::FillDataBlock(const edm::Event      &event, 
                                        const edm::EventSetup &setup)
{
  // Fill the collection.

  //initialize photonfix
  PhotonFix::initialiseGeometry(setup);

  superClusters_->Delete();
  superClusterMap_->Reset();
  superClusterIdMap_->Reset();

  Handle<reco::SuperClusterCollection> hSuperClusterProduct;
  GetProduct(edmName_, hSuperClusterProduct, event);

  Handle<CaloTowerCollection> hCaloTowerProduct;
  GetProduct(caloTowerName_, hCaloTowerProduct, event);

  superClusterMap_->SetEdmProductId(hSuperClusterProduct.id().id());
  const reco::SuperClusterCollection inSuperClusters = *(hSuperClusterProduct.product());  

  // loop through all super clusters
  for (reco::SuperClusterCollection::const_iterator inSC = inSuperClusters.begin(); 
       inSC != inSuperClusters.end(); ++inSC) {

    mithep::SuperCluster *outSC = superClusters_->Allocate();
    new (outSC) mithep::SuperCluster();

    outSC->SetXYZ(inSC->x(),inSC->y(),inSC->z());
    outSC->SetEnergy(inSC->energy());
    outSC->SetRawEnergy(inSC->rawEnergy());
    outSC->SetPreshowerEnergy(inSC->preshowerEnergy());
    outSC->SetPhiWidth(inSC->phiWidth());
    outSC->SetEtaWidth(inSC->etaWidth());

    //Compute Hadronic Energy behind the supercluster (within DR < 0.15)
    EgammaTowerIsolation towerIsoDepth1(0.15,0.,0.,1,hCaloTowerProduct.product()) ;  
    EgammaTowerIsolation towerIsoDepth2(0.15,0.,0.,2,hCaloTowerProduct.product()) ;  
    outSC->SetHcalDepth1Energy(towerIsoDepth1.getTowerESum(&(*inSC)));
    outSC->SetHcalDepth2Energy(towerIsoDepth2.getTowerESum(&(*inSC)));

    // set the seed
    if (basicClusterMap_ && inSC->seed().isNonnull())
      outSC->SetSeed(basicClusterMap_->GetMit(inSC->seed()));
    
    // add basic clusters that belong to this super cluster
    if (basicClusterMap_) {
      for(reco::CaloCluster_iterator bc = inSC->clustersBegin(); bc != inSC->clustersEnd(); ++bc) {
        if (bc->isNonnull())
          outSC->AddCluster(basicClusterMap_->GetMit(*bc));
      }
    }

    //add super cluster det ids to the id map and also fill supercluster-calotower associations
    const std::vector< std::pair<DetId, float> > &pairs = inSC->hitsAndFractions();
    for (std::vector< std::pair<DetId, float> >::const_iterator ipair = pairs.begin();
          ipair < pairs.end(); ++ipair) {

      const DetId &ihit = ipair->first;
    
      if (caloTowerDetIdMap_) {
        if (caloTowerDetIdMap_->HasMit(ihit)) {
          const mithep::CaloTower *matchedTower = caloTowerDetIdMap_->GetMit(ihit);
          if (!outSC->HasTower(matchedTower)) {
            outSC->AddTower(matchedTower);
          }
        }
      }
      
      superClusterIdMap_->Add(ihit,outSC);
    }

    //fill photonfix local coordinates
    PhotonFix pfix(inSC->eta(),inSC->phi());
    
    if (std::abs(inSC->eta())<1.48) {
      outSC->SetEtaC(pfix.etaC());
      outSC->SetEtaS(pfix.etaS());
      outSC->SetEtaM(pfix.etaM());
      outSC->SetPhiC(pfix.phiC());
      outSC->SetPhiS(pfix.phiS());
      outSC->SetPhiM(pfix.phiM());
    }
    else {
      outSC->SetXC(pfix.xC());
      outSC->SetXS(pfix.xS());
      outSC->SetXM(pfix.xM());      
      outSC->SetXZ(pfix.xZ());     
      outSC->SetYC(pfix.yC());
      outSC->SetYS(pfix.yS());
      outSC->SetYM(pfix.yM());      
      outSC->SetYZ(pfix.yZ());         
    }

    // add super cluster to the map
    reco::SuperClusterRef theRef(hSuperClusterProduct, inSC-inSuperClusters.begin());
    superClusterMap_->Add(theRef, outSC);

  }
  superClusters_->Trim();
}
Revision:	1.14
Committed:	Sun Oct 9 23:28:48 2011 UTC (13 years, 6 months ago) by bendavid
Content type:	text/plain
Branch:	MAIN
CVS Tags:	Mit_025e, Mit_025d, Mit_025c, Mit_025b, Mit_025a, Mit_025
Branch point for:	Mit_025c_branch
Changes since 1.13:	+26 -1 lines
Log Message:	store ecal local coordinates
#	User	Rev	Content
1	bendavid	1.14	// $Id: FillerSuperClusters.cc,v 1.13 2010/09/19 23:47:53 bendavid Exp $
2	sixie	1.1
3			#include "MitProd/TreeFiller/interface/FillerSuperClusters.h"
4			#include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
5			#include "DataFormats/EgammaReco/interface/SuperCluster.h"
6			#include "MitAna/DataTree/interface/BasicCluster.h"
7	loizides	1.5	#include "MitAna/DataTree/interface/SuperClusterCol.h"
8	sixie	1.1	#include "MitAna/DataTree/interface/Names.h"
9	loizides	1.5	#include "MitProd/ObjectService/interface/ObjectService.h"
10	sixie	1.10	#include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaTowerIsolation.h"
11	bendavid	1.14	#include "HiggsAnalysis/HiggsToGammaGamma/interface/PhotonFix.h"
12	sixie	1.1
13			using namespace std;
14			using namespace edm;
15			using namespace mithep;
16
17			//--------------------------------------------------------------------------------------------------
18			FillerSuperClusters::FillerSuperClusters(const ParameterSet &cfg, const char *name, bool active) :
19			BaseFiller(cfg,name,active),
20			edmName_(Conf().getUntrackedParameter<string>("edmName","hybridSuperClusters")),
21			mitName_(Conf().getUntrackedParameter<string>("mitName","SuperClusters")),
22			basicClusterMapName_(Conf().getUntrackedParameter<string>("basicClusterMapName",
23			"BasicClusterMap")),
24	bendavid	1.11	caloTowerDetIdMapName_(Conf().getUntrackedParameter<string>("caloTowerDetIdMapName",
25			"CaloTowerDetIdMap")),
26	sixie	1.1	superClusterMapName_(Conf().getUntrackedParameter<string>("superClusterMapName",
27			"SuperClusterMap")),
28	bendavid	1.2	superClusterIdMapName_(Conf().getUntrackedParameter<string>("superClusterIdMapName",
29	loizides	1.7	"SuperClusterIdMap")),
30	sixie	1.10	caloTowerName_(Conf().getUntrackedParameter<string>("caloTowerName","towerMaker")),
31	sixie	1.1	basicClusterMap_(0),
32	bendavid	1.11	caloTowerDetIdMap_(0),
33	sixie	1.1	superClusters_(new mithep::SuperClusterArr(25)),
34	bendavid	1.2	superClusterMap_(new mithep::SuperClusterMap),
35			superClusterIdMap_(new mithep::SuperClusterIdMap)
36	sixie	1.1	{
37			// Constructor.
38			}
39
40			//--------------------------------------------------------------------------------------------------
41			FillerSuperClusters::~FillerSuperClusters()
42			{
43			// Destructor.
44
45			delete superClusters_;
46			delete superClusterMap_;
47	bendavid	1.2	delete superClusterIdMap_;
48	sixie	1.1	}
49
50			//--------------------------------------------------------------------------------------------------
51	bendavid	1.9	void FillerSuperClusters::BookDataBlock(TreeWriter &tws)
52	sixie	1.1	{
53	loizides	1.4	// Add super cluster branch to tree and get pointers to maps.
54
55			tws.AddBranch(mitName_,&superClusters_);
56			OS()->add<SuperClusterArr>(superClusters_,mitName_);
57
58			if (!basicClusterMapName_.empty()) {
59			basicClusterMap_ = OS()->get<BasicClusterMap>(basicClusterMapName_);
60			if (basicClusterMap_)
61			AddBranchDep(mitName_,basicClusterMap_->GetBrName());
62			}
63	sixie	1.1
64	bendavid	1.11	if (!caloTowerDetIdMapName_.empty()) {
65			caloTowerDetIdMap_ = OS()->get<CaloTowerDetIdMap>(caloTowerDetIdMapName_);
66			if (caloTowerDetIdMap_)
67			AddBranchDep(mitName_,caloTowerDetIdMap_->GetBrName());
68			}
69
70	loizides	1.4	if (!superClusterMapName_.empty()) {
71			superClusterMap_->SetBrName(mitName_);
72			OS()->add<SuperClusterMap>(superClusterMap_,superClusterMapName_);
73			}
74			if (!superClusterIdMapName_.empty()) {
75			superClusterIdMap_->SetBrName(mitName_);
76			OS()->add<SuperClusterIdMap>(superClusterIdMap_,superClusterIdMapName_);
77			}
78	sixie	1.1	}
79
80			//--------------------------------------------------------------------------------------------------
81			void FillerSuperClusters::FillDataBlock(const edm::Event &event,
82	loizides	1.4	const edm::EventSetup &setup)
83	sixie	1.1	{
84	loizides	1.4	// Fill the collection.
85
86	bendavid	1.14	//initialize photonfix
87			PhotonFix::initialiseGeometry(setup);
88
89	bendavid	1.3	superClusters_->Delete();
90	sixie	1.1	superClusterMap_->Reset();
91	bendavid	1.2	superClusterIdMap_->Reset();
92	sixie	1.1
93			Handle<reco::SuperClusterCollection> hSuperClusterProduct;
94			GetProduct(edmName_, hSuperClusterProduct, event);
95	sixie	1.10
96			Handle<CaloTowerCollection> hCaloTowerProduct;
97			GetProduct(caloTowerName_, hCaloTowerProduct, event);
98
99	sixie	1.1	superClusterMap_->SetEdmProductId(hSuperClusterProduct.id().id());
100			const reco::SuperClusterCollection inSuperClusters = *(hSuperClusterProduct.product());
101
102			// loop through all super clusters
103			for (reco::SuperClusterCollection::const_iterator inSC = inSuperClusters.begin();
104			inSC != inSuperClusters.end(); ++inSC) {
105
106			mithep::SuperCluster *outSC = superClusters_->Allocate();
107			new (outSC) mithep::SuperCluster();
108
109			outSC->SetXYZ(inSC->x(),inSC->y(),inSC->z());
110			outSC->SetEnergy(inSC->energy());
111			outSC->SetRawEnergy(inSC->rawEnergy());
112			outSC->SetPreshowerEnergy(inSC->preshowerEnergy());
113			outSC->SetPhiWidth(inSC->phiWidth());
114			outSC->SetEtaWidth(inSC->etaWidth());
115
116	sixie	1.10	//Compute Hadronic Energy behind the supercluster (within DR < 0.15)
117			EgammaTowerIsolation towerIsoDepth1(0.15,0.,0.,1,hCaloTowerProduct.product()) ;
118			EgammaTowerIsolation towerIsoDepth2(0.15,0.,0.,2,hCaloTowerProduct.product()) ;
119			outSC->SetHcalDepth1Energy(towerIsoDepth1.getTowerESum(&(*inSC)));
120			outSC->SetHcalDepth2Energy(towerIsoDepth2.getTowerESum(&(*inSC)));
121
122	loizides	1.4	// set the seed
123	sixie	1.1	if (basicClusterMap_ && inSC->seed().isNonnull())
124			outSC->SetSeed(basicClusterMap_->GetMit(inSC->seed()));
125
126	loizides	1.4	// add basic clusters that belong to this super cluster
127	bendavid	1.13	if (basicClusterMap_) {
128			for(reco::CaloCluster_iterator bc = inSC->clustersBegin(); bc != inSC->clustersEnd(); ++bc) {
129			if (bc->isNonnull())
130			outSC->AddCluster(basicClusterMap_->GetMit(*bc));
131			}
132	sixie	1.1	}
133
134	bendavid	1.11	//add super cluster det ids to the id map and also fill supercluster-calotower associations
135	bendavid	1.6	const std::vector< std::pair<DetId, float> > &pairs = inSC->hitsAndFractions();
136			for (std::vector< std::pair<DetId, float> >::const_iterator ipair = pairs.begin();
137			ipair < pairs.end(); ++ipair) {
138	bendavid	1.2
139	bendavid	1.6	const DetId &ihit = ipair->first;
140	bendavid	1.11
141			if (caloTowerDetIdMap_) {
142			if (caloTowerDetIdMap_->HasMit(ihit)) {
143	bendavid	1.12	const mithep::CaloTower *matchedTower = caloTowerDetIdMap_->GetMit(ihit);
144			if (!outSC->HasTower(matchedTower)) {
145			outSC->AddTower(matchedTower);
146			}
147	bendavid	1.11	}
148			}
149
150	bendavid	1.6	superClusterIdMap_->Add(ihit,outSC);
151	bendavid	1.2	}
152	bendavid	1.11
153	bendavid	1.14	//fill photonfix local coordinates
154			PhotonFix pfix(inSC->eta(),inSC->phi());
155
156			if (std::abs(inSC->eta())<1.48) {
157			outSC->SetEtaC(pfix.etaC());
158			outSC->SetEtaS(pfix.etaS());
159			outSC->SetEtaM(pfix.etaM());
160			outSC->SetPhiC(pfix.phiC());
161			outSC->SetPhiS(pfix.phiS());
162			outSC->SetPhiM(pfix.phiM());
163			}
164			else {
165			outSC->SetXC(pfix.xC());
166			outSC->SetXS(pfix.xS());
167			outSC->SetXM(pfix.xM());
168			outSC->SetXZ(pfix.xZ());
169			outSC->SetYC(pfix.yC());
170			outSC->SetYS(pfix.yS());
171			outSC->SetYM(pfix.yM());
172			outSC->SetYZ(pfix.yZ());
173			}
174	bendavid	1.11
175			// add super cluster to the map
176			reco::SuperClusterRef theRef(hSuperClusterProduct, inSC-inSuperClusters.begin());
177			superClusterMap_->Add(theRef, outSC);
178
179	sixie	1.1	}
180			superClusters_->Trim();
181			}