| 1 |
paus |
1.16 |
// $Id: FillerSuperClusters.cc,v 1.15 2012/03/30 01:08:41 paus 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 |
paus |
1.16 |
#include "DataFormats/EgammaReco/interface/PreshowerCluster.h"
|
| 8 |
|
|
#include "DataFormats/EgammaReco/interface/PreshowerClusterFwd.h"
|
| 9 |
loizides |
1.5 |
#include "MitAna/DataTree/interface/SuperClusterCol.h"
|
| 10 |
sixie |
1.1 |
#include "MitAna/DataTree/interface/Names.h"
|
| 11 |
loizides |
1.5 |
#include "MitProd/ObjectService/interface/ObjectService.h"
|
| 12 |
sixie |
1.10 |
#include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaTowerIsolation.h"
|
| 13 |
paus |
1.16 |
#include "RecoEcal/EgammaCoreTools/interface/Mustache.h"
|
| 14 |
|
|
#include "Geometry/CaloEventSetup/interface/CaloTopologyRecord.h"
|
| 15 |
|
|
#include "Geometry/CaloGeometry/interface/CaloGeometry.h"
|
| 16 |
|
|
#include "Geometry/CaloTopology/interface/CaloSubdetectorTopology.h"
|
| 17 |
|
|
#include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
|
| 18 |
|
|
#include "RecoEcal/EgammaCoreTools/interface/EcalClusterLazyTools.h"
|
| 19 |
|
|
#include "Geometry/CaloTopology/interface/EcalPreshowerTopology.h"
|
| 20 |
|
|
#include "Geometry/EcalAlgo/interface/EcalPreshowerGeometry.h"
|
| 21 |
|
|
#include "RecoCaloTools/Navigation/interface/EcalPreshowerNavigator.h"
|
| 22 |
|
|
#include "Geometry/Records/interface/CaloGeometryRecord.h"
|
| 23 |
sixie |
1.1 |
|
| 24 |
|
|
using namespace std;
|
| 25 |
|
|
using namespace edm;
|
| 26 |
|
|
using namespace mithep;
|
| 27 |
|
|
|
| 28 |
|
|
//--------------------------------------------------------------------------------------------------
|
| 29 |
|
|
FillerSuperClusters::FillerSuperClusters(const ParameterSet &cfg, const char *name, bool active) :
|
| 30 |
paus |
1.16 |
BaseFiller (cfg,name,active),
|
| 31 |
|
|
edmName_ (Conf().getUntrackedParameter<string>("edmName","hybridSuperClusters")),
|
| 32 |
|
|
mitName_ (Conf().getUntrackedParameter<string>("mitName","SuperClusters")),
|
| 33 |
|
|
basicClusterMapName_ (Conf().getUntrackedParameter<string>("basicClusterMapName",
|
| 34 |
|
|
"BasicClusterMap")),
|
| 35 |
|
|
psClusterMapName_ (Conf().getUntrackedParameter<string>("psClusterMapName",
|
| 36 |
|
|
"")),
|
| 37 |
bendavid |
1.11 |
caloTowerDetIdMapName_(Conf().getUntrackedParameter<string>("caloTowerDetIdMapName",
|
| 38 |
paus |
1.16 |
"CaloTowerDetIdMap")),
|
| 39 |
|
|
superClusterMapName_ (Conf().getUntrackedParameter<string>("superClusterMapName",
|
| 40 |
|
|
"SuperClusterMap")),
|
| 41 |
bendavid |
1.2 |
superClusterIdMapName_(Conf().getUntrackedParameter<string>("superClusterIdMapName",
|
| 42 |
loizides |
1.7 |
"SuperClusterIdMap")),
|
| 43 |
paus |
1.16 |
caloTowerName_ (Conf().getUntrackedParameter<string>("caloTowerName","towerMaker")),
|
| 44 |
|
|
basicClusterMap_ (0),
|
| 45 |
|
|
psClusterMap_ (0),
|
| 46 |
|
|
caloTowerDetIdMap_ (0),
|
| 47 |
|
|
superClusters_ (new mithep::SuperClusterArr(25)),
|
| 48 |
|
|
superClusterMap_ (new mithep::SuperClusterMap),
|
| 49 |
|
|
superClusterIdMap_ (new mithep::SuperClusterIdMap)
|
| 50 |
sixie |
1.1 |
{
|
| 51 |
|
|
// Constructor.
|
| 52 |
|
|
}
|
| 53 |
|
|
|
| 54 |
|
|
//--------------------------------------------------------------------------------------------------
|
| 55 |
|
|
FillerSuperClusters::~FillerSuperClusters()
|
| 56 |
|
|
{
|
| 57 |
|
|
// Destructor.
|
| 58 |
|
|
|
| 59 |
|
|
delete superClusters_;
|
| 60 |
|
|
delete superClusterMap_;
|
| 61 |
bendavid |
1.2 |
delete superClusterIdMap_;
|
| 62 |
sixie |
1.1 |
}
|
| 63 |
|
|
|
| 64 |
|
|
//--------------------------------------------------------------------------------------------------
|
| 65 |
bendavid |
1.9 |
void FillerSuperClusters::BookDataBlock(TreeWriter &tws)
|
| 66 |
sixie |
1.1 |
{
|
| 67 |
loizides |
1.4 |
// Add super cluster branch to tree and get pointers to maps.
|
| 68 |
|
|
|
| 69 |
|
|
tws.AddBranch(mitName_,&superClusters_);
|
| 70 |
|
|
OS()->add<SuperClusterArr>(superClusters_,mitName_);
|
| 71 |
|
|
|
| 72 |
|
|
if (!basicClusterMapName_.empty()) {
|
| 73 |
|
|
basicClusterMap_ = OS()->get<BasicClusterMap>(basicClusterMapName_);
|
| 74 |
|
|
if (basicClusterMap_)
|
| 75 |
|
|
AddBranchDep(mitName_,basicClusterMap_->GetBrName());
|
| 76 |
|
|
}
|
| 77 |
paus |
1.16 |
|
| 78 |
|
|
if (!psClusterMapName_.empty()) {
|
| 79 |
|
|
psClusterMap_ = OS()->get<PsClusterMap>(psClusterMapName_);
|
| 80 |
|
|
if (psClusterMap_)
|
| 81 |
|
|
AddBranchDep(mitName_,psClusterMap_->GetBrName());
|
| 82 |
|
|
}
|
| 83 |
|
|
|
| 84 |
bendavid |
1.11 |
if (!caloTowerDetIdMapName_.empty()) {
|
| 85 |
|
|
caloTowerDetIdMap_ = OS()->get<CaloTowerDetIdMap>(caloTowerDetIdMapName_);
|
| 86 |
|
|
if (caloTowerDetIdMap_)
|
| 87 |
|
|
AddBranchDep(mitName_,caloTowerDetIdMap_->GetBrName());
|
| 88 |
|
|
}
|
| 89 |
|
|
|
| 90 |
loizides |
1.4 |
if (!superClusterMapName_.empty()) {
|
| 91 |
|
|
superClusterMap_->SetBrName(mitName_);
|
| 92 |
|
|
OS()->add<SuperClusterMap>(superClusterMap_,superClusterMapName_);
|
| 93 |
|
|
}
|
| 94 |
|
|
if (!superClusterIdMapName_.empty()) {
|
| 95 |
|
|
superClusterIdMap_->SetBrName(mitName_);
|
| 96 |
|
|
OS()->add<SuperClusterIdMap>(superClusterIdMap_,superClusterIdMapName_);
|
| 97 |
|
|
}
|
| 98 |
sixie |
1.1 |
}
|
| 99 |
|
|
|
| 100 |
|
|
//--------------------------------------------------------------------------------------------------
|
| 101 |
|
|
void FillerSuperClusters::FillDataBlock(const edm::Event &event,
|
| 102 |
loizides |
1.4 |
const edm::EventSetup &setup)
|
| 103 |
sixie |
1.1 |
{
|
| 104 |
loizides |
1.4 |
// Fill the collection.
|
| 105 |
|
|
|
| 106 |
bendavid |
1.3 |
superClusters_->Delete();
|
| 107 |
sixie |
1.1 |
superClusterMap_->Reset();
|
| 108 |
bendavid |
1.2 |
superClusterIdMap_->Reset();
|
| 109 |
sixie |
1.1 |
|
| 110 |
|
|
Handle<reco::SuperClusterCollection> hSuperClusterProduct;
|
| 111 |
|
|
GetProduct(edmName_, hSuperClusterProduct, event);
|
| 112 |
sixie |
1.10 |
|
| 113 |
|
|
Handle<CaloTowerCollection> hCaloTowerProduct;
|
| 114 |
|
|
GetProduct(caloTowerName_, hCaloTowerProduct, event);
|
| 115 |
|
|
|
| 116 |
sixie |
1.1 |
superClusterMap_->SetEdmProductId(hSuperClusterProduct.id().id());
|
| 117 |
|
|
const reco::SuperClusterCollection inSuperClusters = *(hSuperClusterProduct.product());
|
| 118 |
|
|
|
| 119 |
paus |
1.16 |
//lazy tools for supercluster ecal timing
|
| 120 |
|
|
EcalClusterLazyTools lazyTools(event, setup, edm::InputTag("reducedEcalRecHitsEB"),
|
| 121 |
|
|
edm::InputTag("reducedEcalRecHitsEE"));
|
| 122 |
|
|
|
| 123 |
|
|
//mustache id
|
| 124 |
|
|
reco::Mustache mustache;
|
| 125 |
|
|
|
| 126 |
|
|
//es shape variables
|
| 127 |
|
|
edm::Handle<EcalRecHitCollection> hESRecHits;
|
| 128 |
|
|
event.getByLabel("reducedEcalRecHitsES" , hESRecHits);
|
| 129 |
|
|
|
| 130 |
|
|
edm::ESHandle<CaloGeometry> pGeometry;
|
| 131 |
|
|
setup.get<CaloGeometryRecord>().get(pGeometry);
|
| 132 |
|
|
|
| 133 |
|
|
edm::ESHandle<CaloTopology> pTopology;
|
| 134 |
|
|
setup.get<CaloTopologyRecord>().get(pTopology);
|
| 135 |
|
|
|
| 136 |
|
|
const CaloSubdetectorGeometry *geometryES = pGeometry->getSubdetectorGeometry(DetId::Ecal, EcalPreshower);
|
| 137 |
|
|
EcalPreshowerTopology topology_p(pGeometry);
|
| 138 |
|
|
|
| 139 |
|
|
//map of preshower rechits for shape calculations
|
| 140 |
|
|
std::map<DetId, EcalRecHit> esmap;
|
| 141 |
|
|
if (hESRecHits.isValid()) {
|
| 142 |
|
|
EcalRecHitCollection::const_iterator it;
|
| 143 |
|
|
for (it = hESRecHits->begin(); it != hESRecHits->end(); ++it) {
|
| 144 |
|
|
// remove bad ES rechits
|
| 145 |
|
|
if (it->recoFlag()==1 || it->recoFlag()==14 || (it->recoFlag()<=10 && it->recoFlag()>=5)) continue;
|
| 146 |
|
|
//Make the map of DetID, EcalRecHit pairs
|
| 147 |
|
|
esmap.insert(std::make_pair(it->id(), *it));
|
| 148 |
|
|
}
|
| 149 |
|
|
}
|
| 150 |
|
|
|
| 151 |
sixie |
1.1 |
// loop through all super clusters
|
| 152 |
|
|
for (reco::SuperClusterCollection::const_iterator inSC = inSuperClusters.begin();
|
| 153 |
|
|
inSC != inSuperClusters.end(); ++inSC) {
|
| 154 |
|
|
|
| 155 |
|
|
mithep::SuperCluster *outSC = superClusters_->Allocate();
|
| 156 |
|
|
new (outSC) mithep::SuperCluster();
|
| 157 |
|
|
|
| 158 |
|
|
outSC->SetXYZ(inSC->x(),inSC->y(),inSC->z());
|
| 159 |
|
|
outSC->SetEnergy(inSC->energy());
|
| 160 |
|
|
outSC->SetRawEnergy(inSC->rawEnergy());
|
| 161 |
|
|
outSC->SetPreshowerEnergy(inSC->preshowerEnergy());
|
| 162 |
|
|
outSC->SetPhiWidth(inSC->phiWidth());
|
| 163 |
|
|
outSC->SetEtaWidth(inSC->etaWidth());
|
| 164 |
|
|
|
| 165 |
sixie |
1.10 |
//Compute Hadronic Energy behind the supercluster (within DR < 0.15)
|
| 166 |
|
|
EgammaTowerIsolation towerIsoDepth1(0.15,0.,0.,1,hCaloTowerProduct.product()) ;
|
| 167 |
|
|
EgammaTowerIsolation towerIsoDepth2(0.15,0.,0.,2,hCaloTowerProduct.product()) ;
|
| 168 |
|
|
outSC->SetHcalDepth1Energy(towerIsoDepth1.getTowerESum(&(*inSC)));
|
| 169 |
|
|
outSC->SetHcalDepth2Energy(towerIsoDepth2.getTowerESum(&(*inSC)));
|
| 170 |
|
|
|
| 171 |
paus |
1.16 |
//ecal timing information
|
| 172 |
|
|
outSC->SetTime(lazyTools.SuperClusterTime(*inSC, event));
|
| 173 |
|
|
outSC->SetSeedTime(lazyTools.SuperClusterSeedTime(*inSC));
|
| 174 |
|
|
|
| 175 |
|
|
//preshower shape variables
|
| 176 |
|
|
if (inSC->seed()->hitsAndFractions().at(0).first.subdetId()==EcalEndcap) {
|
| 177 |
|
|
std::vector<float> phoESShape = getESShape(getESHits(inSC->x(), inSC->y(), inSC->z(), esmap, *pGeometry.product(), &topology_p, 0));
|
| 178 |
|
|
|
| 179 |
|
|
if (phoESShape.size()>=2) {
|
| 180 |
|
|
outSC->SetPsEffWidthSigmaXX(phoESShape[0]);
|
| 181 |
|
|
outSC->SetPsEffWidthSigmaYY(phoESShape[1]);
|
| 182 |
|
|
}
|
| 183 |
|
|
|
| 184 |
|
|
}
|
| 185 |
|
|
|
| 186 |
|
|
|
| 187 |
loizides |
1.4 |
// set the seed
|
| 188 |
sixie |
1.1 |
if (basicClusterMap_ && inSC->seed().isNonnull())
|
| 189 |
|
|
outSC->SetSeed(basicClusterMap_->GetMit(inSC->seed()));
|
| 190 |
|
|
|
| 191 |
paus |
1.16 |
// add basic clusters that belong to this super cluster and tag with mustache id
|
| 192 |
bendavid |
1.13 |
if (basicClusterMap_) {
|
| 193 |
paus |
1.16 |
std::vector<reco::CaloCluster> bcs;
|
| 194 |
|
|
std::vector<unsigned int> insidebcs;
|
| 195 |
|
|
std::vector<unsigned int> outsidebcs;
|
| 196 |
bendavid |
1.13 |
for(reco::CaloCluster_iterator bc = inSC->clustersBegin(); bc != inSC->clustersEnd(); ++bc) {
|
| 197 |
paus |
1.16 |
if (bc->isNonnull()) {
|
| 198 |
bendavid |
1.13 |
outSC->AddCluster(basicClusterMap_->GetMit(*bc));
|
| 199 |
paus |
1.16 |
bcs.push_back(**bc);
|
| 200 |
|
|
}
|
| 201 |
|
|
}
|
| 202 |
|
|
|
| 203 |
|
|
mustache.MustacheClust(bcs,insidebcs,outsidebcs);
|
| 204 |
|
|
for (unsigned int iclus=0; iclus<insidebcs.size(); ++iclus) {
|
| 205 |
|
|
const_cast<mithep::BasicCluster*>(outSC->Cluster(insidebcs.at(iclus)))->SetInsideMustache(kTRUE);
|
| 206 |
|
|
}
|
| 207 |
|
|
for (unsigned int iclus=0; iclus<outsidebcs.size(); ++iclus) {
|
| 208 |
|
|
const_cast<mithep::BasicCluster*>(outSC->Cluster(outsidebcs.at(iclus)))->SetInsideMustache(kFALSE);
|
| 209 |
|
|
}
|
| 210 |
|
|
}
|
| 211 |
|
|
|
| 212 |
|
|
// add preshower clusters that belong to this super cluster and tag with mustache id
|
| 213 |
|
|
if (psClusterMap_) {
|
| 214 |
|
|
std::vector<reco::CaloCluster> bcs;
|
| 215 |
|
|
std::vector<unsigned int> insidebcs;
|
| 216 |
|
|
std::vector<unsigned int> outsidebcs;
|
| 217 |
|
|
for(reco::CaloCluster_iterator bc = inSC->preshowerClustersBegin(); bc != inSC->preshowerClustersEnd(); ++bc) {
|
| 218 |
|
|
if (bc->isNonnull()) {
|
| 219 |
|
|
outSC->AddPsClust(psClusterMap_->GetMit(*bc));
|
| 220 |
|
|
bcs.push_back(**bc);
|
| 221 |
|
|
}
|
| 222 |
|
|
}
|
| 223 |
|
|
|
| 224 |
|
|
mustache.MustacheClust(bcs,insidebcs,outsidebcs);
|
| 225 |
|
|
for (unsigned int iclus=0; iclus<insidebcs.size(); ++iclus) {
|
| 226 |
|
|
const_cast<mithep::PsCluster*>(outSC->PsClust(insidebcs.at(iclus)))->SetInsideMustache(kTRUE);
|
| 227 |
|
|
}
|
| 228 |
|
|
for (unsigned int iclus=0; iclus<outsidebcs.size(); ++iclus) {
|
| 229 |
|
|
const_cast<mithep::PsCluster*>(outSC->PsClust(outsidebcs.at(iclus)))->SetInsideMustache(kFALSE);
|
| 230 |
bendavid |
1.13 |
}
|
| 231 |
sixie |
1.1 |
}
|
| 232 |
|
|
|
| 233 |
paus |
1.16 |
//compute preshower energy sums by plane
|
| 234 |
|
|
double psenergyplane1 = 0.;
|
| 235 |
|
|
double psenergyplane2 = 0.;
|
| 236 |
|
|
for(reco::CaloCluster_iterator bc = inSC->preshowerClustersBegin(); bc != inSC->preshowerClustersEnd(); ++bc) {
|
| 237 |
|
|
const reco::PreshowerCluster *pscluster = dynamic_cast<const reco::PreshowerCluster*>(bc->get());
|
| 238 |
|
|
if (!pscluster) continue;
|
| 239 |
|
|
|
| 240 |
|
|
if (pscluster->plane()==1) psenergyplane1 += pscluster->energy();
|
| 241 |
|
|
else if (pscluster->plane()==2) psenergyplane2 += pscluster->energy();
|
| 242 |
|
|
}
|
| 243 |
|
|
outSC->SetPreshowerEnergyPlane1(psenergyplane1);
|
| 244 |
|
|
outSC->SetPreshowerEnergyPlane2(psenergyplane2);
|
| 245 |
|
|
|
| 246 |
bendavid |
1.11 |
//add super cluster det ids to the id map and also fill supercluster-calotower associations
|
| 247 |
bendavid |
1.6 |
const std::vector< std::pair<DetId, float> > &pairs = inSC->hitsAndFractions();
|
| 248 |
|
|
for (std::vector< std::pair<DetId, float> >::const_iterator ipair = pairs.begin();
|
| 249 |
|
|
ipair < pairs.end(); ++ipair) {
|
| 250 |
bendavid |
1.2 |
|
| 251 |
bendavid |
1.6 |
const DetId &ihit = ipair->first;
|
| 252 |
bendavid |
1.11 |
|
| 253 |
|
|
if (caloTowerDetIdMap_) {
|
| 254 |
|
|
if (caloTowerDetIdMap_->HasMit(ihit)) {
|
| 255 |
bendavid |
1.12 |
const mithep::CaloTower *matchedTower = caloTowerDetIdMap_->GetMit(ihit);
|
| 256 |
|
|
if (!outSC->HasTower(matchedTower)) {
|
| 257 |
|
|
outSC->AddTower(matchedTower);
|
| 258 |
|
|
}
|
| 259 |
bendavid |
1.11 |
}
|
| 260 |
|
|
}
|
| 261 |
|
|
|
| 262 |
bendavid |
1.6 |
superClusterIdMap_->Add(ihit,outSC);
|
| 263 |
bendavid |
1.2 |
}
|
| 264 |
bendavid |
1.11 |
|
| 265 |
|
|
// add super cluster to the map
|
| 266 |
|
|
reco::SuperClusterRef theRef(hSuperClusterProduct, inSC-inSuperClusters.begin());
|
| 267 |
|
|
superClusterMap_->Add(theRef, outSC);
|
| 268 |
|
|
|
| 269 |
sixie |
1.1 |
}
|
| 270 |
|
|
superClusters_->Trim();
|
| 271 |
|
|
}
|
| 272 |
paus |
1.16 |
|
| 273 |
|
|
//horrible code below copied from globe for preshower cluster shape calculations
|
| 274 |
|
|
std::vector<float> FillerSuperClusters::getESHits(double X, double Y, double Z, std::map<DetId, EcalRecHit> rechits_map, const CaloGeometry& geometry, CaloSubdetectorTopology *topology_p, int row) {
|
| 275 |
|
|
std::vector<float> esHits;
|
| 276 |
|
|
|
| 277 |
|
|
const GlobalPoint point(X,Y,Z);
|
| 278 |
|
|
|
| 279 |
|
|
const CaloSubdetectorGeometry *geometry_p ;
|
| 280 |
|
|
geometry_p = geometry.getSubdetectorGeometry (DetId::Ecal,EcalPreshower) ;
|
| 281 |
|
|
|
| 282 |
|
|
DetId esId1 = (dynamic_cast<const EcalPreshowerGeometry*>(geometry_p))->getClosestCellInPlane(point, 1);
|
| 283 |
|
|
DetId esId2 = (dynamic_cast<const EcalPreshowerGeometry*>(geometry_p))->getClosestCellInPlane(point, 2);
|
| 284 |
|
|
ESDetId esDetId1 = (esId1 == DetId(0)) ? ESDetId(0) : ESDetId(esId1);
|
| 285 |
|
|
ESDetId esDetId2 = (esId2 == DetId(0)) ? ESDetId(0) : ESDetId(esId2);
|
| 286 |
|
|
|
| 287 |
|
|
std::map<DetId, EcalRecHit>::iterator it;
|
| 288 |
|
|
ESDetId next;
|
| 289 |
|
|
ESDetId strip1;
|
| 290 |
|
|
ESDetId strip2;
|
| 291 |
|
|
|
| 292 |
|
|
strip1 = esDetId1;
|
| 293 |
|
|
strip2 = esDetId2;
|
| 294 |
|
|
|
| 295 |
|
|
EcalPreshowerNavigator theESNav1(strip1, topology_p);
|
| 296 |
|
|
theESNav1.setHome(strip1);
|
| 297 |
|
|
|
| 298 |
|
|
EcalPreshowerNavigator theESNav2(strip2, topology_p);
|
| 299 |
|
|
theESNav2.setHome(strip2);
|
| 300 |
|
|
|
| 301 |
|
|
if (row == 1) {
|
| 302 |
|
|
if (strip1 != ESDetId(0)) strip1 = theESNav1.north();
|
| 303 |
|
|
if (strip2 != ESDetId(0)) strip2 = theESNav2.east();
|
| 304 |
|
|
} else if (row == -1) {
|
| 305 |
|
|
if (strip1 != ESDetId(0)) strip1 = theESNav1.south();
|
| 306 |
|
|
if (strip2 != ESDetId(0)) strip2 = theESNav2.west();
|
| 307 |
|
|
}
|
| 308 |
|
|
|
| 309 |
|
|
// Plane 1
|
| 310 |
|
|
if (strip1 == ESDetId(0)) {
|
| 311 |
|
|
for (int i=0; i<31; ++i) esHits.push_back(0);
|
| 312 |
|
|
} else {
|
| 313 |
|
|
|
| 314 |
|
|
it = rechits_map.find(strip1);
|
| 315 |
|
|
if (it->second.energy() > 1.0e-10 && it != rechits_map.end()) esHits.push_back(it->second.energy());
|
| 316 |
|
|
else esHits.push_back(0);
|
| 317 |
|
|
//cout<<"center : "<<strip1<<" "<<it->second.energy()<<endl;
|
| 318 |
|
|
|
| 319 |
|
|
// east road
|
| 320 |
|
|
for (int i=0; i<15; ++i) {
|
| 321 |
|
|
next = theESNav1.east();
|
| 322 |
|
|
if (next != ESDetId(0)) {
|
| 323 |
|
|
it = rechits_map.find(next);
|
| 324 |
|
|
if (it->second.energy() > 1.0e-10 && it != rechits_map.end()) esHits.push_back(it->second.energy());
|
| 325 |
|
|
else esHits.push_back(0);
|
| 326 |
|
|
//cout<<"east "<<i<<" : "<<next<<" "<<it->second.energy()<<endl;
|
| 327 |
|
|
} else {
|
| 328 |
|
|
for (int j=i; j<15; j++) esHits.push_back(0);
|
| 329 |
|
|
break;
|
| 330 |
|
|
//cout<<"east "<<i<<" : "<<next<<" "<<0<<endl;
|
| 331 |
|
|
}
|
| 332 |
|
|
}
|
| 333 |
|
|
|
| 334 |
|
|
// west road
|
| 335 |
|
|
theESNav1.setHome(strip1);
|
| 336 |
|
|
theESNav1.home();
|
| 337 |
|
|
for (int i=0; i<15; ++i) {
|
| 338 |
|
|
next = theESNav1.west();
|
| 339 |
|
|
if (next != ESDetId(0)) {
|
| 340 |
|
|
it = rechits_map.find(next);
|
| 341 |
|
|
if (it->second.energy() > 1.0e-10 && it != rechits_map.end()) esHits.push_back(it->second.energy());
|
| 342 |
|
|
else esHits.push_back(0);
|
| 343 |
|
|
//cout<<"west "<<i<<" : "<<next<<" "<<it->second.energy()<<endl;
|
| 344 |
|
|
} else {
|
| 345 |
|
|
for (int j=i; j<15; j++) esHits.push_back(0);
|
| 346 |
|
|
break;
|
| 347 |
|
|
//cout<<"west "<<i<<" : "<<next<<" "<<0<<endl;
|
| 348 |
|
|
}
|
| 349 |
|
|
}
|
| 350 |
|
|
}
|
| 351 |
|
|
|
| 352 |
|
|
if (strip2 == ESDetId(0)) {
|
| 353 |
|
|
for (int i=0; i<31; ++i) esHits.push_back(0);
|
| 354 |
|
|
} else {
|
| 355 |
|
|
|
| 356 |
|
|
it = rechits_map.find(strip2);
|
| 357 |
|
|
if (it->second.energy() > 1.0e-10 && it != rechits_map.end()) esHits.push_back(it->second.energy());
|
| 358 |
|
|
else esHits.push_back(0);
|
| 359 |
|
|
//cout<<"center : "<<strip2<<" "<<it->second.energy()<<endl;
|
| 360 |
|
|
|
| 361 |
|
|
// north road
|
| 362 |
|
|
for (int i=0; i<15; ++i) {
|
| 363 |
|
|
next = theESNav2.north();
|
| 364 |
|
|
if (next != ESDetId(0)) {
|
| 365 |
|
|
it = rechits_map.find(next);
|
| 366 |
|
|
if (it->second.energy() > 1.0e-10 && it != rechits_map.end()) esHits.push_back(it->second.energy());
|
| 367 |
|
|
else esHits.push_back(0);
|
| 368 |
|
|
//cout<<"north "<<i<<" : "<<next<<" "<<it->second.energy()<<endl;
|
| 369 |
|
|
} else {
|
| 370 |
|
|
for (int j=i; j<15; j++) esHits.push_back(0);
|
| 371 |
|
|
break;
|
| 372 |
|
|
//cout<<"north "<<i<<" : "<<next<<" "<<0<<endl;
|
| 373 |
|
|
}
|
| 374 |
|
|
}
|
| 375 |
|
|
|
| 376 |
|
|
// south road
|
| 377 |
|
|
theESNav2.setHome(strip2);
|
| 378 |
|
|
theESNav2.home();
|
| 379 |
|
|
for (int i=0; i<15; ++i) {
|
| 380 |
|
|
next = theESNav2.south();
|
| 381 |
|
|
if (next != ESDetId(0)) {
|
| 382 |
|
|
it = rechits_map.find(next);
|
| 383 |
|
|
if (it->second.energy() > 1.0e-10 && it != rechits_map.end()) esHits.push_back(it->second.energy());
|
| 384 |
|
|
else esHits.push_back(0);
|
| 385 |
|
|
//cout<<"south "<<i<<" : "<<next<<" "<<it->second.energy()<<endl;
|
| 386 |
|
|
} else {
|
| 387 |
|
|
for (int j=i; j<15; j++) esHits.push_back(0);
|
| 388 |
|
|
break;
|
| 389 |
|
|
//cout<<"south "<<i<<" : "<<next<<" "<<0<<endl;
|
| 390 |
|
|
}
|
| 391 |
|
|
}
|
| 392 |
|
|
}
|
| 393 |
|
|
|
| 394 |
|
|
return esHits;
|
| 395 |
|
|
}
|
| 396 |
|
|
|
| 397 |
|
|
std::vector<float> FillerSuperClusters::getESShape(std::vector<float> ESHits0)
|
| 398 |
|
|
{
|
| 399 |
|
|
std::vector<float> esShape;
|
| 400 |
|
|
|
| 401 |
|
|
const int nBIN = 21;
|
| 402 |
|
|
float esRH_F[nBIN];
|
| 403 |
|
|
float esRH_R[nBIN];
|
| 404 |
|
|
for (int idx=0; idx<nBIN; idx++) {
|
| 405 |
|
|
esRH_F[idx] = 0.;
|
| 406 |
|
|
esRH_R[idx] = 0.;
|
| 407 |
|
|
}
|
| 408 |
|
|
|
| 409 |
|
|
for(int ibin=0; ibin<((nBIN+1)/2); ibin++) {
|
| 410 |
|
|
if (ibin==0) {
|
| 411 |
|
|
esRH_F[(nBIN-1)/2] = ESHits0[ibin];
|
| 412 |
|
|
esRH_R[(nBIN-1)/2] = ESHits0[ibin+31];
|
| 413 |
|
|
} else {
|
| 414 |
|
|
esRH_F[(nBIN-1)/2+ibin] = ESHits0[ibin];
|
| 415 |
|
|
esRH_F[(nBIN-1)/2-ibin] = ESHits0[ibin+15];
|
| 416 |
|
|
esRH_R[(nBIN-1)/2+ibin] = ESHits0[ibin+31];
|
| 417 |
|
|
esRH_R[(nBIN-1)/2-ibin] = ESHits0[ibin+31+15];
|
| 418 |
|
|
}
|
| 419 |
|
|
}
|
| 420 |
|
|
|
| 421 |
|
|
// ---- Effective Energy Deposit Width ---- //
|
| 422 |
|
|
double EffWidthSigmaXX = 0.;
|
| 423 |
|
|
double EffWidthSigmaYY = 0.;
|
| 424 |
|
|
double totalEnergyXX = 0.;
|
| 425 |
|
|
double totalEnergyYY = 0.;
|
| 426 |
|
|
double EffStatsXX = 0.;
|
| 427 |
|
|
double EffStatsYY = 0.;
|
| 428 |
|
|
for (int id_X=0; id_X<21; id_X++) {
|
| 429 |
|
|
totalEnergyXX += esRH_F[id_X];
|
| 430 |
|
|
EffStatsXX += esRH_F[id_X]*(id_X-10)*(id_X-10);
|
| 431 |
|
|
totalEnergyYY += esRH_R[id_X];
|
| 432 |
|
|
EffStatsYY += esRH_R[id_X]*(id_X-10)*(id_X-10);
|
| 433 |
|
|
}
|
| 434 |
|
|
EffWidthSigmaXX = (totalEnergyXX>0.) ? sqrt(fabs(EffStatsXX / totalEnergyXX)) : 0.;
|
| 435 |
|
|
EffWidthSigmaYY = (totalEnergyYY>0.) ? sqrt(fabs(EffStatsYY / totalEnergyYY)) : 0.;
|
| 436 |
|
|
|
| 437 |
|
|
esShape.push_back(EffWidthSigmaXX);
|
| 438 |
|
|
esShape.push_back(EffWidthSigmaYY);
|
| 439 |
|
|
|
| 440 |
|
|
return esShape;
|
| 441 |
|
|
}
|
| 442 |
|
|
|
