Skip to content

Commit

Permalink
add RecHits validator CPU/GPU
Browse files Browse the repository at this point in the history
  • Loading branch information
Bruno Alves committed Aug 28, 2020
1 parent 03f23b3 commit 9bca0cb
Show file tree
Hide file tree
Showing 3 changed files with 272 additions and 0 deletions.
83 changes: 83 additions & 0 deletions UserCode/RecHitsValidator/plugins/Validator.cc
View file
@@ -0,0 +1,83 @@
#include "RecoLocalCalo/HGCalRecProducers/plugins/HeterogeneousHGCalHEFRecHitProducer.h"

#include "DetectorDescription/OfflineDBLoader/interface/GeometryInfoDump.h"
#include "Geometry/Records/interface/IdealGeometryRecord.h"
#include "Geometry/HGCalGeometry/interface/HGCalGeometry.h"
#include "SimG4CMS/Calo/interface/CaloHitID.h"

#include "DetectorDescription/Core/interface/DDFilter.h"
#include "DetectorDescription/Core/interface/DDFilteredView.h"
#include "DetectorDescription/Core/interface/DDSolid.h"

#include "DataFormats/GeometryVector/interface/Basic3DVector.h"

#include "CLHEP/Units/GlobalSystemOfUnits.h"

#include "FWCore/ServiceRegistry/interface/Service.h"
#include "FWCore/Utilities/interface/Exception.h"
#include "CommonTools/UtilAlgos/interface/TFileService.h"

HGCalMaskResolutionAna::HGCalMaskResolutionAna( const edm::ParameterSet &ps ) :
mc_( consumes<edm::HepMCProduct>(edm::InputTag("generatorSmeared")) ),
recHitsTokens_( {consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>("recHitsCEEToken")),
consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>("recHitsHSiToken")),
consumes<HGCRecHitCollection>(ps.getParameter<edm::InputTag>("recHitsHScToken"))} ),
genParticles_( consumes<std::vector<reco::GenParticle>>(edm::InputTag("genParticles"))),
treename_("tree")
{

edm::Service<TFileService> fs;
for(size_t idet = 0; idet<nSubDets_; ++idet) {
zhist.push_back(fs->make<TH1F>(("z"+treename_).c_str(), ("z"+treename_).c_str(), 2000, -10, 650));
xhist.push_back(fs->make<TH1F>(("x"+treename_).c_str(), ("x"+treename_).c_str(), 1000, -40, 40));
yhist.push_back(fs->make<TH1F>(("y"+treename_).c_str(), ("y"+treename_).c_str(), 1000, -40, 40));
zsidehist.push_back(fs->make<TH1F>(("zside"+treename_).c_str(), ("zside"+treename_).c_str(), 4, -1.1, 1.1));
ehist.push_back(fs->make<TH1F>(("e"+treename_).c_str(), ("e"+treename_).c_str(), 2000, -10, 650));
layerhist.push_back(fs->make<TH1F>(("layer"+treename_).c_str(), ("layer"+treename_).c_str(), 102, 0, 52));
offsetlayerhist.push_back(fs->make<TH1F>(("offsetlayer"+treename_).c_str(), ("offsetlayer"+treename_).c_str(), 102, 0, 52));
}

tree_ = fs->make<TTree>((treenames_[0]+"_"+treenames_[1]+"_"+treenames_[2]).c_str(),
(treenames_[0]+"_"+treenames_[1]+"_"+treenames_[2]).c_str());
tree_->Branch("Hits", "std::vector<SlimmedHit>", &slimmedHits_);
tree_->Branch("ROIs", "std::vector<SlimmedROI>", &slimmedROIs_);
tree_->Branch("GenVertex", "TLorentzVector", &genVertex_);
}

HGCalMaskResolutionAna::~HGCalMaskResolutionAna()
{
}

void HGCalMaskResolutionAna::endJob()
{
}

void HGCalMaskResolutionAna::analyze(const edm::Event &iEvent, const edm::EventSetup &iSetup)
{
recHitTools_.getEventSetup(iSetup);

std::vector< edm::Handle<HGCRecHitCollection> > recHitsHandles(nSubDets_);
std::vector< edm::ESHandle<HGCalGeometry> > geomHandle(nSubDets_);

//read the generated primary vertex
edm::Handle<edm::HepMCProduct> mcHandle;
iEvent.getByToken(mc_, mcHandle);
HepMC::GenVertex *primaryVertex = *(mcHandle)->GetEvent()->vertices_begin();

///CEE, HFE (Silicon), HBF (Scintillator)///
for(size_t idet=0; idet<nSubDets_; ++idet) {
iSetup.get<IdealGeometryRecord>().get(geometrySource_[idet], geomHandle[idet]);

//collect rec hits in regions of interest
iEvent.getByToken(recHitsTokens_[idet], recHitsHandles[idet]);
int nMatched(0), nUnmatched(0);
for(size_t i=0; i<recHitsHandles[idet]->size(); i++) {
const HGCRecHit &h = recHitsHandles[idet]->operator[](i);
const HGCalDetId did = h.detid();
}

tree_->Fill();
}

//define this as a plug-in
DEFINE_FWK_MODULE(HGCalMaskResolutionAna);
64 changes: 64 additions & 0 deletions UserCode/RecHitsValidator/plugins/Validator.h
View file
@@ -0,0 +1,64 @@
#ifndef _HGCalMaskResolutionAna_h_
#define _HGCalMaskResolutionAna_h_

#include "FWCore/Framework/interface/Frameworkfwd.h"
#include "FWCore/Framework/interface/EDAnalyzer.h"

#include "FWCore/Framework/interface/Event.h"
#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/Framework/interface/ESTransientHandle.h"
#include "FWCore/Framework/interface/MakerMacros.h"
#include "FWCore/ServiceRegistry/interface/Service.h"

#include "CommonTools/UtilAlgos/interface/TFileService.h"

#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"

#include "DataFormats/GeometryVector/interface/GlobalPoint.h"
#include "DataFormats/HepMCCandidate/interface/GenParticle.h"
#include "DataFormats/HepMCCandidate/interface/GenParticleFwd.h"
#include "DataFormats/HGCRecHit/interface/HGCRecHit.h"
#include "DataFormats/HGCRecHit/interface/HGCRecHitCollections.h"

#include "RecoLocalCalo/HGCalRecAlgos/interface/RecHitTools.h"

#include "SimDataFormats/Track/interface/SimTrack.h"
#include "SimDataFormats/Track/interface/SimTrackContainer.h"
#include "SimDataFormats/Vertex/interface/SimVertex.h"
#include "SimDataFormats/Vertex/interface/SimVertexContainer.h"

#include "TTree.h"
#include "TH1F.h"

#include <iostream>
#include <string>

class HeterogeneousHGCalRecHitsValidator : public edm::EDAnalyzer
{
public:
explicit HeterogeneousHGCalRecHitsValidator( const edm::ParameterSet& );
~HeterogeneousHGCalRecHitsValidator();
virtual void analyze( const edm::Event&, const edm::EventSetup& );
void endJob();

private:
edm::EDGetTokenT<edm::HepMCProduct> mc_;
std::vector< edm::EDGetTokenT<HGCRecHitCollection> > recHitsTokens_;
edm::EDGetTokenT<reco::GenParticleCollection> genParticles_;

hgcal::RecHitTools recHitTools_;

std::vector< TH1F* > zhist;
std::vector< TH1F* > xhist;
std::vector< TH1F* > yhist;
std::vector< TH1F* > zsidehist;
std::vector< TH1F* > ehist;
std::vector< TH1F* > layerhist;
std::vector< TH1F* > offsetlayerhist;

TTree* tree_;
std::string treename_;
};


#endif
125 changes: 125 additions & 0 deletions UserCode/RecHitsValidator/python/Validator_cfg.py
View file
@@ -0,0 +1,125 @@
import os, sys, glob
import FWCore.ParameterSet.Config as cms
from Configuration.StandardSequences.Eras import eras

#arguments parsing
from FWCore.ParameterSet.VarParsing import VarParsing
F = VarParsing('analysis')
F.register('pu',
0,
F.multiplicity.singleton,
F.varType.bool,
"Whether to run with pile-up.")
F.register('fidx',
0,
F.multiplicity.singleton,
F.varType.int,
"Which file index to consider.")
F.register('mask',
-1,
F.multiplicity.singleton,
F.varType.int,
"Mask to be used. Accepted values: 3, 4, 5 or 6.")
F.register('samples',
'',
F.multiplicity.singleton,
F.varType.string,
'Which samples to use. Inner ("inner"), outer ("outer"), or both ("all").')
F.parseArguments()
print("********************")
print("Input arguments:")
for k,v in F.__dict__["_singletons"].items():
print("{}: {}".format(k,v))
print("********************")

#package loading
process = cms.Process("postRECO", eras.Phase2C8)
process.load('Configuration.StandardSequences.Services_cff')
process.load("FWCore.MessageLogger.MessageLogger_cfi")
process.load('SimGeneral.HepPDTESSource.pythiapdt_cfi')
process.load('Configuration.StandardSequences.MagneticField_cff')
process.load('Configuration.EventContent.EventContent_cff')
process.load('Configuration.StandardSequences.FrontierConditions_GlobalTag_cff')
process.load('Configuration.Geometry.GeometryExtended2023D28_cff')
process.load('Configuration.Geometry.GeometryExtended2023D28Reco_cff')
"""
process.MessageLogger = cms.Service("MessageLogger",
destinations = cms.untracked.vstring(
'detailedInfo'),
detailedInfo = cms.untracked.PSet(
threshold = cms.untracked.string('INFO'),
default = cms.untracked.PSet(
limit = cms.untracked.int32(-1))),
debugModules = cms.untracked.vstring('*'))
"""

from Configuration.AlCa.GlobalTag import GlobalTag
process.GlobalTag = GlobalTag(process.GlobalTag, 'auto:phase2_realistic', '')
from RecoLocalCalo.HGCalRecProducers.HGCalRecHit_cfi import dEdX_weights_v10

process.maxEvents = cms.untracked.PSet( input = cms.untracked.int32(-1) )

indir1 = "/eos/cms/store/cmst3/group/hgcal/CMG_studies/Production/FlatRandomEGunProducer_PionGun_SciStudies_bfontana_inner_20190716/RECO/"
indir2 = "/eos/cms/store/cmst3/group/hgcal/CMG_studies/Production/FlatRandomEGunProducer_PionGun_SciStudies_bfontana_outer_20190716/RECO/"
indir3 = "/eos/cms/store/cmst3/group/hgcal/CMG_studies/Production/FlatRandomEGunProducer_PionGun_CalibrationStudies_bfontana_central_20190822/RECO/"

glob1 = glob.glob(os.path.join(indir1,"*.root"))
glob2 = glob.glob(os.path.join(indir2,"*.root"))
glob3 = glob.glob(os.path.join(indir3,"*.root"))
if F.samples == 'all':
glob_tot = glob1 + glob2
elif F.samples == 'inner':
glob_tot = glob1
elif F.samples == 'outer':
glob_tot = glob2
elif F.samples == 'central':
glob_tot = glob3
else:
raise ValueError('Insert a valid "samples" option!')
print("Total number of files: ", len(glob_tot))
fNames = ["file:" + it for it in glob_tot][F.fidx]
print("this file: ", fNames)

if isinstance(fNames,list):
process.source = cms.Source("PoolSource",
fileNames = cms.untracked.vstring(fNames),
duplicateCheckMode = cms.untracked.string("noDuplicateCheck"))
else:
process.source = cms.Source("PoolSource",
fileNames = cms.untracked.vstring(fNames),
duplicateCheckMode = cms.untracked.string("noDuplicateCheck"))

process.RecHitsMasked = cms.EDProducer('HGCalMaskProd',
recHitsCEEToken = cms.InputTag('HGCalRecHit', 'HGCEERecHits'),
recHitsHSiToken = cms.InputTag('HGCalRecHit', 'HGCHEFRecHits'),
Mask = cms.uint32(F.mask))

process.an = cms.EDAnalyzer("HGCalMaskResolutionAna",
recHitsCEEToken = cms.InputTag('HGCalRecHit','HGCEERecHits'),
recHitsHSiToken = cms.InputTag('HGCalRecHit','HGCHEFRecHits'),
recHitsHScToken = cms.InputTag('HGCalRecHit','HGCHEBRecHits'),
thicknessCorrection = cms.vdouble(0.781, 0.776, 0.769),
dEdXWeights = dEdX_weights_v10,
geometrySource = cms.vstring('HGCalEESensitive',
'HGCalHESiliconSensitive',
'HGCalHEScintillatorSensitive'),
distancesSR1 = cms.vdouble(30., 30., 30.),
distancesSR2 = cms.vdouble(40., 40., 40.),
distancesSR3 = cms.vdouble(50., 50., 50.),
nControlRegions = cms.int32(2), #5 for photons
particle = cms.string('pion'),
byClosest = cms.bool(False))

process.an_mask = process.an.clone(
recHitsCEEToken = cms.InputTag("RecHitsMasked", "HGCEERecHits"),
recHitsHSiToken = cms.InputTag("RecHitsMasked", "HGCHSiRecHits"))

pu_str = "pu" if F.pu else "nopu"
fileName = str(F.fidx)+"_mask"+str(F.mask)+"_"+F.samples+"_"+pu_str
process.TFileService = cms.Service("TFileService",
fileName = cms.string(fileName+".root"))
#fileName = fileName + '_out'
#process.out = cms.OutputModule("PoolOutputModule",
#fileName = cms.untracked.string(fileName+".root"))
process.p = cms.Path(process.RecHitsMasked * process.an * process.an_mask)
#process.outpath = cms.EndPath(process.out)

0 comments on commit 9bca0cb

Please sign in to comment.