RECO Step working for HGCal + Basic CMS PFlow Hooks

DataFormats/CaloRecHit/interface/CaloID.h

@@ -29,6 +29,9 @@ namespace reco {
       DET_HF_EM,
       DET_HF_HAD,
       DET_HO,
+      DET_HGC_ECAL,
+      DET_HGC_HCALB,
+      DET_HGC_HCALF,
       DET_NONE
     };
 

DataFormats/ParticleFlowReco/interface/PFBlockElement.h

@@ -43,7 +43,9 @@ namespace reco {
       HFHAD=9,
       SC=10,
       HO=11,
-      kNBETypes=12
+      ECAL_HGC=12,
+      HCAL_HGC=13,
+      kNBETypes=14
     };
 
     enum TrackType {

DataFormats/ParticleFlowReco/interface/PFLayer.h

@@ -37,7 +37,10 @@ class PFLayer {
               HCAL_BARREL2 = 2,
               HCAL_ENDCAP  = 3,
               HF_EM        = 11,
-              HF_HAD       = 12};
+              HF_HAD       = 12,
+              HGC_ECAL     = 13,
+              HGC_HCALF    = 14,
+	      HGC_HCALB    = 15  };
 
   static reco::CaloID  toCaloID( Layer layer);
 

DataFormats/ParticleFlowReco/src/PFLayer.cc

@@ -20,6 +20,9 @@ CaloID  PFLayer::toCaloID( Layer layer) {
   case HCAL_ENDCAP   : return  CaloID(CaloID::DET_HCAL_ENDCAP); 
   case HF_EM         : return  CaloID(CaloID::DET_HF_EM); 
   case HF_HAD        : return  CaloID(CaloID::DET_HF_HAD); 
+  case HGC_ECAL      : return  CaloID(CaloID::DET_HGC_ECAL);
+  case HGC_HCALF     : return  CaloID(CaloID::DET_HGC_HCALF);
+  case HGC_HCALB     : return  CaloID(CaloID::DET_HGC_HCALB);
   default            : return  CaloID();
   }
 }
@@ -41,7 +44,10 @@ PFLayer::Layer   PFLayer::fromCaloID( const CaloID& id) {
   case CaloID::DET_HCAL_ENDCAP   : return  HCAL_ENDCAP;
   case CaloID::DET_HF_EM         : return  HF_EM;
   case CaloID::DET_HF_HAD        : return  HF_HAD;
-  case CaloID::DET_HO            : return  HCAL_BARREL2; 
+  case CaloID::DET_HO            : return  HCAL_BARREL2;
+  case CaloID::DET_HGC_ECAL      : return  HGC_ECAL;
+  case CaloID::DET_HGC_HCALF     : return  HGC_HCALF;
+  case CaloID::DET_HGC_HCALB     : return  HGC_HCALB;
   default                        : return  NONE;
   }
 }

RecoEgamma/EgammaIsolationAlgos/plugins/EgammaRecHitExtractor.cc

@@ -124,7 +124,7 @@ reco::IsoDeposit EgammaRecHitExtractor::deposit(const edm::Event & iEvent,
   const CaloSubdetectorGeometry* endcapgeom = caloGeom->getSubdetectorGeometry(DetId::Ecal,EcalEndcap);
   if( endcapgeom==nullptr ) endcapgeom=caloGeom->getSubdetectorGeometry(DetId::Ecal,EcalShashlik); //Shervin
 
-  assert(endcapgeom!=nullptr);
+  //assert(endcapgeom!=nullptr);
 
   static std::string metname = "EgammaIsolationAlgos|EgammaRecHitExtractor";
 
@@ -160,7 +160,7 @@ reco::IsoDeposit EgammaRecHitExtractor::deposit(const edm::Event & iEvent,
     collect(deposit, sc, barrelgeom, caloGeom, *barrelEcalRecHitsH, sevLevel, true);
   } 
 
-  if ((!inBarrel) || tryBoth_) {
+  if ( endcapgeom != nullptr && ( (!inBarrel) || tryBoth_ ) ) {
     collect(deposit, sc, endcapgeom, caloGeom, *endcapEcalRecHitsH, sevLevel, false);
   }
 

RecoEgamma/EgammaIsolationAlgos/src/EgammaRecHitIsolation.cc

@@ -61,7 +61,7 @@ EgammaRecHitIsolation::EgammaRecHitIsolation (double extRadius,
     subdet_[0] = caloGeom->getSubdetectorGeometry(DetId::Ecal,EcalBarrel);
     subdet_[1] = caloGeom->getSubdetectorGeometry(DetId::Ecal,EcalEndcap);
     if(subdet_[1]==NULL) subdet_[1] = caloGeom->getSubdetectorGeometry(DetId::Ecal,EcalShashlik);
-    assert(subdet_[1]!=NULL);
+    //assert(subdet_[1]!=NULL);
 }
 
 EgammaRecHitIsolation::~EgammaRecHitIsolation ()
@@ -84,6 +84,7 @@ double EgammaRecHitIsolation::getSum_(const reco::Candidate* emObject,bool retur
     std::vector< std::pair<DetId, float> >::const_iterator rhIt;
 
     for(int subdetnr=0; subdetnr<=1 ; subdetnr++){  // look in barrel and endcap
+      if( !subdet_[subdetnr] ) continue;
       CaloSubdetectorGeometry::DetIdSet chosen = subdet_[subdetnr]->getCells(pclu,extRadius_);// select cells around cluster
       CaloRecHitMetaCollectionV::const_iterator j = caloHits_->end();
       for (CaloSubdetectorGeometry::DetIdSet::const_iterator  i = chosen.begin ();i != chosen.end (); ++i){ //loop selected cells
@@ -199,6 +200,7 @@ double EgammaRecHitIsolation::getSum_(const reco::SuperCluster* sc, bool returnE
 
 
     for(int subdetnr=0; subdetnr<=1 ; subdetnr++){  // look in barrel and endcap
+      if( !subdet_[subdetnr] ) continue;
       CaloSubdetectorGeometry::DetIdSet chosen = subdet_[subdetnr]->getCells(pclu,extRadius_);// select cells around cluster
       CaloRecHitMetaCollectionV::const_iterator j=caloHits_->end();
       for (CaloSubdetectorGeometry::DetIdSet::const_iterator  i = chosen.begin ();i!= chosen.end ();++i){//loop selected cells

RecoLocalCalo/CaloTowersCreator/src/CaloTowersCreationAlgo.cc

@@ -1304,13 +1304,15 @@ GlobalPoint CaloTowersCreationAlgo::hadShwPosFromCells(DetId frontCellId, DetId
 
     GlobalPoint point = frontCellGeometry->getPosition();
 
-    CaloCellGeometry::CornersVec cv = backCellGeometry->getCorners();
-
-    GlobalPoint backPoint = GlobalPoint(0.25 * (cv[4].x() + cv[5].x() + cv[6].x() + cv[7].x()),
-      0.25 * (cv[4].y() + cv[5].y() + cv[6].y() + cv[7].y()),
-      0.25 * (cv[4].z() + cv[5].z() + cv[6].z() + cv[7].z()));
+    if( backCellGeometry ) {
+      CaloCellGeometry::CornersVec cv = backCellGeometry->getCorners();
+
+      GlobalPoint backPoint = GlobalPoint(0.25 * (cv[4].x() + cv[5].x() + cv[6].x() + cv[7].x()),
+					  0.25 * (cv[4].y() + cv[5].y() + cv[6].y() + cv[7].y()),
+					  0.25 * (cv[4].z() + cv[5].z() + cv[6].z() + cv[7].z()));
 
-    point += fracDepth * (backPoint - point);
+      point += fracDepth * (backPoint - point);
+    }
 
     return point;
 }

RecoLocalCalo/Configuration/python/RecoLocalCalo_cff.py

@@ -17,10 +17,13 @@
 #HCAL reconstruction
 from RecoLocalCalo.Configuration.hcalLocalReco_cff import *
 from RecoLocalCalo.Configuration.hcalGlobalReco_cff import *
+# 
+# HGCAL reconstruction
+from RecoLocalCalo.Configuration.hgcalLocalRecoSequence_cff import *
 #
 # sequence CaloLocalReco and CaloGlobalReco
 #
-calolocalreco = cms.Sequence(ecalLocalRecoSequence+hcalLocalRecoSequence)
+calolocalreco = cms.Sequence(ecalLocalRecoSequence+hcalLocalRecoSequence+hgcalLocalRecoSequence)
 caloglobalreco = cms.Sequence(hcalGlobalRecoSequence)
 
 from RecoLocalCalo.HcalRecProducers.HcalHitSelection_cfi import *

RecoLocalCalo/Configuration/python/hgcalLocalReco_EventContent_cff.py

@@ -4,13 +4,17 @@
 #Full Event content 
 hgcalLocalRecoFEVT = cms.PSet(
     outputCommands = cms.untracked.vstring(
-        'keep *_HGCalUncalibRecHit_*_*', 
-        'keep *_HGCRecHit_*_*')
+        'keep *_*_HGCalUncalibRecHit_*',
+        'keep *_*_HGCHEBRecHits_*',
+        'keep *_*_HGCHEFRecHits_*',        
+        'keep *_*_HGCEERecHits_*')
 )
 #RECO content
 hgcalLocalRecoRECO = cms.PSet(
     outputCommands = cms.untracked.vstring(
-        'keep *_HGCRecHit_*_*',
+        'keep *_*_HGCHEBRecHits_*',
+        'keep *_*_HGCHEFRecHits_*',
+        'keep *_*_HGCEERecHits_*'        
         )
 )
 #AOD content

RecoLocalCalo/HGCalRecAlgos/interface/HGCalRecHitSimpleAlgo.h

@@ -17,8 +17,8 @@ class HGCalRecHitSimpleAlgo : public HGCalRecHitAbsAlgo {
  public:
   // default ctor
   HGCalRecHitSimpleAlgo() {
-    adcToGeVConstant_ = -1;
-    adcToGeVConstantIsSet_ = false;
+    adcToGeVConstant_ = 1;
+    adcToGeVConstantIsSet_ = true;
   }
 
   virtual void setADCToGeVConstant(const float& value) {

RecoLocalCalo/HGCalRecProducers/plugins/BuildFile.xml

@@ -1,5 +1,6 @@
 <use   name="Geometry/CaloGeometry"/>
 <use   name="RecoLocalCalo/HGCalRecProducers"/>
+<use   name="RecoLocalCalo/HGCalRecAlgos"/>
 <use   name="FWCore/Framework"/>
 <use   name="FWCore/ParameterSet"/>
 <use   name="DataFormats/HGCDigi"/>

RecoLocalCalo/HGCalRecProducers/python/HGCalUncalibRecHit_cfi.py

@@ -2,11 +2,11 @@
 
 # HGCAL producer of rechits starting from digis
 HGCalUncalibRecHit = cms.EDProducer("HGCalUncalibRecHitProducer",
-    HGCEEdigiCollection = cms.InputTag("HGCalDigis","HGCEEDigis"),
+    HGCEEdigiCollection = cms.InputTag("mix","HGCDigisEE"),
     HGCEEhitCollection = cms.string('HGCEEUncalibRecHits'),
-    HGCHEFdigiCollection = cms.InputTag("HGCalDigis","HGCHEFDigis"),
+    HGCHEFdigiCollection = cms.InputTag("mix","HGCDigisHEfront"),
     HGCHEFhitCollection = cms.string('HGCHEFUncalibRecHits'),
-    HGCHEBdigiCollection = cms.InputTag("HGCalDigis","HGCHEBDigis"),
+    HGCHEBdigiCollection = cms.InputTag("mix","HGCDigisHEback"),
     HGCHEBhitCollection = cms.string('HGCHEBUncalibRecHits'),                                    
     algo = cms.string("HGCalUncalibRecHitWorkerWeights")
 )

RecoParticleFlow/Configuration/python/RecoParticleFlow_EventContent_cff.py

@@ -10,7 +10,7 @@
     outputCommands = cms.untracked.vstring(
     'drop CaloTowersSorted_towerMakerPF_*_*',
     #'keep recoPFRecHits_*_Cleaned_*',
-    'keep recoPFRecHits_particleFlowClusterECAL_Cleaned_*',
+    'keep recoPFRecHits_particleFlowClusterECAL_Cleaned_*',    
     'keep recoPFRecHits_particleFlowClusterHCAL_Cleaned_*',
     'keep recoPFRecHits_particleFlowClusterHO_Cleaned_*',
     'keep recoPFRecHits_particleFlowClusterHFEM_Cleaned_*',
@@ -20,8 +20,11 @@
     'keep recoPFRecHits_particleFlowRecHitHO_Cleaned_*',
     'keep recoPFRecHits_particleFlowRecHitHCAL_Cleaned_*',
     'keep recoPFRecHits_particleFlowRecHitPS_Cleaned_*',
+    'keep recoPFRecHits_particleFlowRecHitEK_*_*',
+    'keep recoPFRecHits_particleFlowRecHitHGC*_*_*',
     #'keep recoPFClusters_*_*_*',
-    'keep recoPFClusters_particleFlowClusterECAL_*_*',    
+    'keep recoPFClusters_particleFlowClusterECAL_*_*',
+    'keep recoPFClusters_particleFlowClusterHGC*_*_*',
     'keep recoPFClusters_particleFlowClusterHCAL_*_*',
     'keep recoPFClusters_particleFlowClusterHO_*_*',
     'keep recoPFClusters_particleFlowClusterHFEM_*_*',
@@ -62,13 +65,16 @@
     'keep recoPFRecHits_particleFlowRecHitHCAL_Cleaned_*',
     'keep recoPFRecHits_particleFlowRecHitHO_Cleaned_*',
     'keep recoPFRecHits_particleFlowRecHitPS_Cleaned_*',
+    'keep recoPFRecHits_particleFlowRecHitEK_Cleaned_*',
+    'keep recoPFRecHits_particleFlowRecHitHGC*_Cleaned_*',
     #'keep recoPFClusters_*_*_*',
     'keep recoPFClusters_particleFlowClusterECAL_*_*',
     'keep recoPFClusters_particleFlowClusterHCAL_*_*',
-    'keep recoPFClusters_particleFlowClusterHO_*_*',
+    'keep recoPFClusters_particleFlowClusterHO_*_*',    
     #'keep recoPFClusters_particleFlowClusterHFEM_*_*',
     #'keep recoPFClusters_particleFlowClusterHFHAD_*_*',
     'keep recoPFClusters_particleFlowClusterPS_*_*',
+    'keep recoPFClusters_particleFlowClusterHGC*_*_*',
     #'keep recoPFBlocks_*_*_*',
     'keep recoPFBlocks_particleFlowBlock_*_*',
     #'keep recoPFCandidates_*_*_*',
@@ -107,6 +113,8 @@
     'keep recoPFRecHits_particleFlowRecHitHCAL_Cleaned_*',
     'keep recoPFRecHits_particleFlowRecHitHO_Cleaned_*',
     'keep recoPFRecHits_particleFlowRecHitPS_Cleaned_*',
+    'keep recoPFRecHits_particleFlowRecHitEK_Cleaned_*',
+    'keep recoPFRecHits_particleFlowRecHitHGC*_Cleaned_*',
     #'keep recoPFCandidates_*_*_*',
     'keep recoCaloClusters_particleFlowEGamma_*_*',
     'keep recoSuperClusters_particleFlowEGamma_*_*',    

RecoParticleFlow/PFClusterProducer/interface/PFHGCalRecHitCreator.h

@@ -0,0 +1,128 @@
+#ifndef RecoParticleFlow_PFClusterProducer_PFHGCalRecHitCreator_h
+#define RecoParticleFlow_PFClusterProducer_PFHGCalRecHitCreator_h
+
+#include "RecoParticleFlow/PFClusterProducer/interface/PFRecHitCreatorBase.h"
+
+#include "Geometry/FCalGeometry/interface/HGCalGeometry.h"
+#include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"
+
+#include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
+#include "Geometry/CaloGeometry/interface/CaloGeometry.h"
+#include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
+#include "Geometry/CaloGeometry/interface/TruncatedPyramid.h"
+#include "Geometry/Records/interface/CaloGeometryRecord.h"
+
+#include "Geometry/EcalAlgo/interface/EcalEndcapGeometry.h"
+#include "Geometry/EcalAlgo/interface/EcalBarrelGeometry.h"
+#include "Geometry/CaloTopology/interface/EcalEndcapTopology.h"
+#include "Geometry/CaloTopology/interface/EcalBarrelTopology.h"
+#include "Geometry/CaloTopology/interface/EcalPreshowerTopology.h"
+#include "RecoCaloTools/Navigation/interface/CaloNavigator.h"
+
+template <typename DET,PFLayer::Layer Layer,ForwardSubdetector subdet>
+  class PFHGCalRecHitCreator :  public  PFRecHitCreatorBase {
+
+ public:  
+  PFHGCalRecHitCreator(const edm::ParameterSet& iConfig,edm::ConsumesCollector& iC):
+    PFRecHitCreatorBase(iConfig,iC)
+    {
+      recHitToken_ = iC.consumes<HGCRecHitCollection>(iConfig.getParameter<edm::InputTag>("src"));
+      geometryInstance_ = iConfig.getParameter<std::string>("geometryInstance");
+    }
+
+    void importRecHits(std::auto_ptr<reco::PFRecHitCollection>&out,std::auto_ptr<reco::PFRecHitCollection>& cleaned ,const edm::Event& iEvent,const edm::EventSetup& iSetup) {
+
+      for (unsigned int i=0;i<qualityTests_.size();++i) {
+	qualityTests_.at(i)->beginEvent(iEvent,iSetup);
+      }
+
+      edm::Handle<HGCRecHitCollection> recHitHandle;
+      iEvent.getByToken(recHitToken_,recHitHandle);
+      const HGCRecHitCollection& rechits = *recHitHandle;
+
+      edm::ESHandle<HGCalGeometry> geoHandle;
+      iSetup.get<IdealGeometryRecord>().get(geometryInstance_,geoHandle);
+      const HGCalGeometry& hgcGeo = *geoHandle;
+
+      unsigned skipped_rechits = 0;
+      for (unsigned int i=0;i<rechits.size();++i) {
+	const HGCRecHit& hgrh = rechits[i];
+	const DET detid(hgrh.detid());
+
+	if( subdet != detid.subdet() ) {
+	  std::cout << "subdet expected: " << subdet 
+		    << " subdet gotten: " << detid.subdet() << std::endl;
+	}
+	assert(subdet == detid.subdet() && "Incorrect subdetector for this importer!");
+	double energy = hgrh.energy();
+	double time = hgrh.time();	
+
+	const FlatTrd *thisCell = 
+	  static_cast<const FlatTrd*>(hgcGeo.getGeometry(detid));	
+
+	// find rechit geometry
+	if(!thisCell) {
+	  LogDebug("PFHGCalRecHitCreator")
+	    <<"warning detid "<<detid.rawId()
+	    <<" not found in geometry"<<std::endl;
+	  ++skipped_rechits;
+	  continue;
+	}
+
+	const auto& position = thisCell->getPosition();
+	//std::cout << "geometry cell position: " << position << std::endl;
+
+	reco::PFRecHit rh( detid.rawId(),Layer,
+			   energy, 
+			   position.x(), position.y(), position.z(), 
+			   0, 0, 0 ); 
+
+	const CaloCellGeometry::CornersVec& corners = thisCell->getCorners();
+	assert( corners.size() == 8 );
+
+	rh.setNECorner( corners[0].x(), corners[0].y(),  corners[0].z() );
+	rh.setSECorner( corners[1].x(), corners[1].y(),  corners[1].z() );
+	rh.setSWCorner( corners[2].x(), corners[2].y(),  corners[2].z() );
+	rh.setNWCorner( corners[3].x(), corners[3].y(),  corners[3].z() );
+
+	bool rcleaned = false;
+	bool keep=true;
+
+	//Apply Q tests
+	for (unsigned int i=0;i<qualityTests_.size();++i) {
+	  if (!qualityTests_.at(i)->test(rh,hgrh,rcleaned)) {
+	    keep = false;	    
+	  }
+	}
+
+	if(keep) {
+	  rh.setTime(time);
+	  out->push_back(rh);
+	}
+	else if (rcleaned) 
+	  cleaned->push_back(rh);
+      }
+      LogDebug("HGCalRecHitCreator") 
+	<<  "Skipped " << skipped_rechits 
+	<< " out of " << rechits.size() << " rechits!" << std::endl;
+      edm::LogInfo("HGCalRecHitCreator")
+	<< "Created " << out->size() << " PFRecHits!" << std::endl;
+    }
+
+
+
+ protected:
+  edm::EDGetTokenT<HGCRecHitCollection> recHitToken_;
+  std::string geometryInstance_;
+
+};
+
+#include "DataFormats/ForwardDetId/interface/HGCEEDetId.h"
+#include "DataFormats/ForwardDetId/interface/HGCHEDetId.h"
+
+typedef PFHGCalRecHitCreator<HGCEEDetId,PFLayer::HGC_ECAL,HGCEE> PFHGCEERecHitCreator;
+typedef PFHGCalRecHitCreator<HGCHEDetId,PFLayer::HGC_HCALF,HGCHEF> PFHGCHEFRecHitCreator;
+typedef PFHGCalRecHitCreator<HGCHEDetId,PFLayer::HGC_HCALB,HGCHEB> PFHGCHEBRecHitCreator;
+
+
+#endif

RecoParticleFlow/PFClusterProducer/interface/PFRecHitQTestBase.h

@@ -21,6 +21,7 @@
 #include "DataFormats/HcalRecHit/interface/HFRecHit.h"
 #include "DataFormats/HcalRecHit/interface/HORecHit.h"
 #include "DataFormats/HcalRecHit/interface/HBHERecHit.h"
+#include "DataFormats/HGCRecHit/interface/HGCRecHit.h"
 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
 #include "CondFormats/HcalObjects/interface/HcalChannelQuality.h"
 #include "CondFormats/DataRecord/interface/HcalChannelQualityRcd.h"
@@ -41,6 +42,7 @@ class PFRecHitQTestBase {
   virtual bool test( reco::PFRecHit& ,const HFRecHit&,bool&)=0;
   virtual bool test( reco::PFRecHit& ,const HORecHit&,bool&)=0;
   virtual bool test( reco::PFRecHit& ,const CaloTower&,bool&)=0;
+  virtual bool test( reco::PFRecHit& ,const HGCRecHit&,bool&)=0;
 };