Skip to content

Commit

Permalink
Merge pull request cms-sw#25 from bendavid/WmassNanoProd_10_6_26_trac…
Browse files Browse the repository at this point in the history
…krefit_v709

Update track fits for nano production (v709)
  • Loading branch information
kdlong authored Dec 16, 2022
2 parents 39fc5c5 + 0b4e7b6 commit 9a6c375
Show file tree
Hide file tree
Showing 27 changed files with 8,616 additions and 21,854 deletions.
3,282 changes: 0 additions & 3,282 deletions Analysis/HitAnalyzer/plugins/ResidualGlobalCorrectionMaker.cc

This file was deleted.

6,077 changes: 2,270 additions & 3,807 deletions Analysis/HitAnalyzer/plugins/ResidualGlobalCorrectionMakerBase.cc

Large diffs are not rendered by default.

254 changes: 111 additions & 143 deletions Analysis/HitAnalyzer/plugins/ResidualGlobalCorrectionMakerBase.h

Large diffs are not rendered by default.

4,099 changes: 1,529 additions & 2,570 deletions Analysis/HitAnalyzer/plugins/ResidualGlobalCorrectionMakerG4e.cc

Large diffs are not rendered by default.

2,430 changes: 0 additions & 2,430 deletions Analysis/HitAnalyzer/plugins/ResidualGlobalCorrectionMakerSim.cc

This file was deleted.

2,459 changes: 0 additions & 2,459 deletions Analysis/HitAnalyzer/plugins/ResidualGlobalCorrectionMakerSimG4e.cc

This file was deleted.

2,320 changes: 0 additions & 2,320 deletions Analysis/HitAnalyzer/plugins/ResidualGlobalCorrectionMakerTwoTrack.cc

This file was deleted.

1,650 changes: 685 additions & 965 deletions Analysis/HitAnalyzer/plugins/ResidualGlobalCorrectionMakerTwoTrackG4e.cc

Large diffs are not rendered by default.

Original file line number Diff line number Diff line change
Expand Up @@ -10,6 +10,7 @@
#include "OAEParametrizedMagneticField.h"
#include "ParabolicParametrizedMagneticField.h"
#include "PolyFit2DParametrizedMagneticField.h"
#include "PolyFit3DParametrizedMagneticField.h"

#include "FWCore/Utilities/interface/Exception.h"

Expand All @@ -32,7 +33,8 @@ ParametrizedMagneticFieldFactory::get(string version, const ParameterSet& parame
return result;
} else if (version=="PolyFit3D") {
// V. Maroussov polynomial fit to mapping data
throw cms::Exception("InvalidParameter")<<"PolyFit3D is not supported anymore";
std::auto_ptr<MagneticField> result( new PolyFit3DParametrizedMagneticField(parameters));
return result;
} else if (version=="Parabolic"){
// FIXME implement configurable parameters to be passed to ctor
// vector<double> params = parameters.getParameter<vdouble>("parameters");
Expand Down
Original file line number Diff line number Diff line change
@@ -0,0 +1,92 @@
/** \file
*
* \author N. Amapane
*/

#include "PolyFit3DParametrizedMagneticField.h"
#include <FWCore/ParameterSet/interface/ParameterSet.h>
#include <FWCore/MessageLogger/interface/MessageLogger.h>

#include "BFit3D.h"


using namespace std;
using namespace magfieldparam;

PolyFit3DParametrizedMagneticField::PolyFit3DParametrizedMagneticField(double bVal) :
theParam(new BFit3D())
{
theParam->SetField(bVal);
}


PolyFit3DParametrizedMagneticField::PolyFit3DParametrizedMagneticField(const edm::ParameterSet& parameters) : theParam(new BFit3D()) {
theParam->SetField(parameters.getParameter<double>("BValue"));

// Additional options (documentation by Vassili):

// By default, the package accepts signed value of "r". That means,
// one can cross r=0 and orientation of the coordinate "orts"
// e_r and e_phi will not be flipped over.
// In other words for an r<0 the e_r points inward, in the direction r=0.
// This is a "natural" mode. However, the default behavior may be
// changed by the call

// theParam->UseSignedRad(false);

// In that case with crossing of r=0 e_r and e_phi will be flipped in
// such a way that e_r always points outward. In other words instead of
// (r<0, phi) the (abs(r), phi+PI) will be used in this mode.

// The expansion coefficients for a nominal field in between measurement
// field values (2.0T, 3.5T, 3.8T and 4.0T) by default are calculated by
// means of a linear piecewise interpolation. Another provided
// interpolation mode is cubic spline. This mode can be switched
// on by the call:

// theParam->UseSpline(true);

// From practical point of view the use of spline interpolation doesn't
// change much, but it makes the coefficients' behavior a bit more
// physical at very low/high field values.
}


PolyFit3DParametrizedMagneticField::~PolyFit3DParametrizedMagneticField() {
delete theParam;
}


GlobalVector
PolyFit3DParametrizedMagneticField::inTesla(const GlobalPoint& gp) const {

if (isDefined(gp)) {
return inTeslaUnchecked(gp);
} else {
edm::LogWarning("MagneticField|FieldOutsideValidity") << " Point " << gp << " is outside the validity region of PolyFit3DParametrizedMagneticField";
return GlobalVector();
}
}

GlobalVector
PolyFit3DParametrizedMagneticField::inTeslaUnchecked(const GlobalPoint& gp) const {
double Br, Bz, Bphi;
theParam->GetField(gp.perp()/100., gp.z()/100., gp.phi(),
Br, Bz, Bphi);

double cosphi = cos(gp.phi());
double sinphi = sin(gp.phi());

return GlobalVector(Br*cosphi - Bphi*sinphi,
Br*sinphi + Bphi*cosphi,
Bz);
}

bool
PolyFit3DParametrizedMagneticField::isDefined(const GlobalPoint& gp) const {
double z = fabs(gp.z());
double r = gp.perp();
//"rectangle" |z|<3.5, r<1.9 _except_ the "corners" |z|+2.5*r>6.7, everything in meters
if (z>350. || r>190 || z+2.5*r>670.) return false;
return true;
}
Original file line number Diff line number Diff line change
@@ -0,0 +1,40 @@
#ifndef PolyFit3DParametrizedMagneticField_h
#define PolyFit3DParametrizedMagneticField_h

/** \class PolyFit3DParametrizedMagneticField
*
* Magnetic Field engine wrapper for V. Maroussov's 3D parametrization
* of the MT data.
*
* \author N. Amapane
*/

#include "MagneticField/Engine/interface/MagneticField.h"

namespace edm { class ParameterSet; }
namespace magfieldparam { class BFit3D; }


class PolyFit3DParametrizedMagneticField : public MagneticField {
public:
/// Constructor. Fitted bVal for the nominal currents are:
/// 2.0216; 3.5162; 3.8114; 4.01242188708911
PolyFit3DParametrizedMagneticField(double bVal = 3.8114);

/// Constructor. Parameters taken from a PSet
PolyFit3DParametrizedMagneticField(const edm::ParameterSet& parameters);

/// Destructor
virtual ~PolyFit3DParametrizedMagneticField();

GlobalVector inTesla (const GlobalPoint& gp) const;

GlobalVector inTeslaUnchecked (const GlobalPoint& gp) const;

bool isDefined(const GlobalPoint& gp) const;

private:
magfieldparam::BFit3D* theParam;
};
#endif

Original file line number Diff line number Diff line change
Expand Up @@ -84,10 +84,10 @@ class FlattenedValueMapVectorTableProducer : public edm::stream::EDProducer<> {
iEvent.getByToken(intVecMaps_[i], vmap);
const auto& results = readVals(*vmap, objs, sizes);

auto intsizetab = std::make_unique<nanoaod::FlatTable>(objs.size(), this->name_+intNames_[i]+"_Counts", false, false);
auto intsizetab = std::make_unique<nanoaod::FlatTable>(objs.size(), this->name_ + "_" + intNames_[i]+"_Counts", false, false);
intsizetab->template addColumn<int>("", sizes, "Number of entries per object", nanoaod::FlatTable::IntColumn, -1);
intsizetab->setDoc(doc_);
auto intvectab = std::make_unique<nanoaod::FlatTable>(results.size(), this->name_+intNames_[i], false, false);
auto intvectab = std::make_unique<nanoaod::FlatTable>(results.size(), this->name_+ "_" + intNames_[i], false, false);
intvectab->template addColumn<int>("Vals", results, intDocs_[i], nanoaod::FlatTable::IntColumn, intPrecisions_[i]);
intvectab->setDoc(doc_);

Expand All @@ -100,10 +100,10 @@ class FlattenedValueMapVectorTableProducer : public edm::stream::EDProducer<> {
iEvent.getByToken(floatVecMaps_[i], vmap);
const auto& results = readVals(*vmap, objs, sizes);

auto floatsizetab = std::make_unique<nanoaod::FlatTable>(objs.size(), this->name_+floatNames_[i]+"_Counts", false, false);
auto floatsizetab = std::make_unique<nanoaod::FlatTable>(objs.size(), this->name_ + "_" + floatNames_[i]+"_Counts", false, false);
floatsizetab->template addColumn<int>("", sizes, "Number of entries per object", nanoaod::FlatTable::IntColumn, -1);
floatsizetab->setDoc(doc_);
auto floatvectab = std::make_unique<nanoaod::FlatTable>(results.size(), this->name_+floatNames_[i], false, false);
auto floatvectab = std::make_unique<nanoaod::FlatTable>(results.size(), this->name_ + "_" + floatNames_[i], false, false);
floatvectab->template addColumn<float>("Vals", results, floatDocs_[i], nanoaod::FlatTable::FloatColumn, floatPrecisions_[i]);
floatvectab->setDoc(doc_);

Expand Down
80 changes: 58 additions & 22 deletions PhysicsTools/NanoAOD/python/muons_cff.py
Original file line number Diff line number Diff line change
Expand Up @@ -128,36 +128,48 @@
trackrefit = cms.EDProducer('ResidualGlobalCorrectionMakerG4e',
src = cms.InputTag("tracksfrommuons"),
fitFromGenParms = cms.bool(False),
fitFromSimParms = cms.bool(False),
fillTrackTree = cms.bool(False),
fillGrads = cms.bool(False),
fillJac = cms.bool(False),
fillRunTree = cms.bool(False),
doGen = cms.bool(False),
doSim = cms.bool(False),
requireGen = cms.bool(False),
doMuons = cms.bool(False),
doMuonAssoc = cms.bool(True),
doTrigger = cms.bool(False),
doRes = cms.bool(False),
bsConstraint = cms.bool(False),
applyHitQuality = cms.bool(True),
corFile = cms.string(""),
useIdealGeometry = cms.bool(False),
corFiles = cms.vstring(),
)

trackrefitbs = cms.EDProducer('ResidualGlobalCorrectionMakerG4e',
trackrefitideal = cms.EDProducer('ResidualGlobalCorrectionMakerG4e',
src = cms.InputTag("tracksfrommuons"),
fitFromGenParms = cms.bool(False),
fitFromSimParms = cms.bool(False),
fillTrackTree = cms.bool(False),
fillGrads = cms.bool(False),
fillJac = cms.bool(False),
fillRunTree = cms.bool(False),
doGen = cms.bool(False),
doSim = cms.bool(False),
requireGen = cms.bool(False),
doMuons = cms.bool(False),
doMuonAssoc = cms.bool(True),
bsConstraint = cms.bool(True),
doTrigger = cms.bool(False),
doRes = cms.bool(False),
bsConstraint = cms.bool(False),
applyHitQuality = cms.bool(True),
corFile = cms.string(""),
useIdealGeometry = cms.bool(True),
corFiles = cms.vstring(),
)

mergedGlobalIdxs = cms.EDProducer("GlobalIdxProducer",
src0 = cms.InputTag("trackrefit", "globalIdxs"),
src1 = cms.InputTag("trackrefitbs", "globalIdxs")
src1 = cms.InputTag("trackrefitideal", "globalIdxs")
)

muonTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
Expand Down Expand Up @@ -229,37 +241,57 @@
cvhPhi = ExtVar(cms.InputTag("trackrefit:corPhi"), float, doc="Refitted track phi", precision=12),
cvhCharge = ExtVar(cms.InputTag("trackrefit:corCharge"), int, doc="Refitted track charge"),
cvhEdmval = ExtVar(cms.InputTag("trackrefit:edmval"), float, doc="Refitted estimated distance to minimum", precision=10),
cvhbsPt = ExtVar(cms.InputTag("trackrefitbs:corPt"), float, doc="Refitted track pt (with bs constraint)", precision=-1),
cvhbsEta = ExtVar(cms.InputTag("trackrefitbs:corEta"), float, doc="Refitted track eta (with bs constraint)", precision=12),
cvhbsPhi = ExtVar(cms.InputTag("trackrefitbs:corPhi"), float, doc="Refitted track phi (with bs constraint)", precision=12),
cvhbsCharge = ExtVar(cms.InputTag("trackrefitbs:corCharge"), int, doc="Refitted track charge (with bs constraint)"),
cvhbsEdmval = ExtVar(cms.InputTag("trackrefitbs:edmval"), float, doc="Refitted estimated distance to minimum (with bs constraint)", precision=10),
),
)

muonIdealTable = cms.EDProducer("SimpleCandidateFlatTableProducer",
src = muonTable.src,
cut = muonTable.cut,
name = muonTable.name,
singleton = cms.bool(False), # the number of entries is variable
extension = cms.bool(True),
variables = cms.PSet(),
externalVariables = cms.PSet(
cvhidealPt = ExtVar(cms.InputTag("trackrefitideal:corPt"), float, doc="Refitted track pt (with ideal geometry)", precision=-1),
cvhidealEta = ExtVar(cms.InputTag("trackrefitideal:corEta"), float, doc="Refitted track eta (with ideal geometry)", precision=12),
cvhidealPhi = ExtVar(cms.InputTag("trackrefitideal:corPhi"), float, doc="Refitted track phi (with ideal geometry)", precision=12),
cvhidealCharge = ExtVar(cms.InputTag("trackrefitideal:corCharge"), int, doc="Refitted track charge (with ideal geometry)"),
cvhidealEdmval = ExtVar(cms.InputTag("trackrefitideal:edmval"), float, doc="Refitted estimated distance to minimum (with ideal geometry)", precision=10),
),
)

muonExternalVecVarsTable = cms.EDProducer("FlattenedCandValueMapVectorTableProducer",
name = cms.string(muonTable.name.value()+"_cvh"),
name = cms.string(muonTable.name.value()),
src = muonTable.src,
cut = muonTable.cut,
doc = muonTable.doc,
variables = cms.PSet(
# can you declare a max number of bits here? technically for the moment this needs 16 bits, but might eventually need 17 or 18
mergedGlobalIdxs = ExtVar(cms.InputTag("mergedGlobalIdxs"), "std::vector<int>", doc="Indices for correction parameters"),
cvhmergedGlobalIdxs = ExtVar(cms.InputTag("trackrefit:globalIdxs"), "std::vector<int>", doc="Indices for correction parameters"),
# optimal precision tbd, but presumably can work the same way as for scalar floats
JacRef = ExtVar(cms.InputTag("trackrefit:jacRef"), "std::vector<float>", doc="jacobian for corrections", precision = 12),
MomCov = ExtVar(cms.InputTag("trackrefit:momCov"), "std::vector<float>", doc="covariance matrix for qop, lambda, phi", precision = 12),
cvhJacRef = ExtVar(cms.InputTag("trackrefit:jacRef"), "std::vector<float>", doc="jacobian for corrections", precision = 12),
cvhMomCov = ExtVar(cms.InputTag("trackrefit:momCov"), "std::vector<float>", doc="covariance matrix for qop, lambda, phi", precision = 12),
)
)

muonIdealExternalVecVarsTable = cms.EDProducer("FlattenedCandValueMapVectorTableProducer",
name = cms.string(muonTable.name.value()),
src = muonTable.src,
cut = muonTable.cut,
doc = muonTable.doc,
variables = cms.PSet(
# can you declare a max number of bits here? technically for the moment this needs 16 bits, but might eventually need 17 or 18
cvhmergedGlobalIdxs = ExtVar(cms.InputTag("mergedGlobalIdxs"), "std::vector<int>", doc="Indices for correction parameters"),
# optimal precision tbd, but presumably can work the same way as for scalar floats
bsJacRef = ExtVar(cms.InputTag("trackrefitbs:jacRef"), "std::vector<float>", doc="Jacobian for corrections (with bs constraint)", precision = 12),
bsMomCov = ExtVar(cms.InputTag("trackrefitbs:momCov"), "std::vector<float>", doc="covariance matrix for qop, lambda, phi (with bs constraint)", precision = 12),
cvhidealJacRef = ExtVar(cms.InputTag("trackrefitideal:jacRef"), "std::vector<float>", doc="jacobian for corrections (with ideal geometry)", precision = 12),
cvhidealMomCov = ExtVar(cms.InputTag("trackrefitideal:momCov"), "std::vector<float>", doc="covariance matrix for qop, lambda, phi (with ideal geometry)", precision = 12),
)
)


for modifier in run2_miniAOD_80XLegacy, run2_nanoAOD_94X2016, run2_nanoAOD_94XMiniAODv1, run2_nanoAOD_94XMiniAODv2, run2_nanoAOD_LowPU:
modifier.toModify(muonTable.variables, puppiIsoId = None, softMva = None)

run2_nanoAOD_LowPU.toModify(muonTable, externalVariables = cms.PSet())

run2_nanoAOD_102Xv1.toModify(muonTable.variables, puppiIsoId = None)

(run2_nanoAOD_106Xv1 & ~run2_nanoAOD_devel).toModify(muonTable.variables, isStandalone = None)
Expand All @@ -285,9 +317,13 @@
docString = cms.string("MC matching to status==1 muons"),
)

muonSequence = cms.Sequence(slimmedMuonsUpdated+isoForMu + ptRatioRelForMu + slimmedMuonsWithUserData + finalMuons + finalLooseMuons )
muonMC = cms.Sequence(muonsMCMatchForTable + muonMCTable)
muonTables = cms.Sequence(muonFSRphotons + muonFSRassociation + muonMVATTH + muonMVALowPt + geopro + tracksfrommuons + trackrefit + trackrefitbs + mergedGlobalIdxs + muonTable + muonExternalVecVarsTable + fsrTable)
muonSequence = cms.Sequence(slimmedMuonsUpdated+isoForMu + ptRatioRelForMu + slimmedMuonsWithUserData + finalMuons + finalLooseMuons)
muonMC = cms.Sequence(trackrefitideal + mergedGlobalIdxs + muonsMCMatchForTable + muonMCTable + muonIdealTable + muonIdealExternalVecVarsTable)
muonTables = cms.Sequence(muonFSRphotons + muonFSRassociation + muonMVATTH + muonMVALowPt + geopro + tracksfrommuons + trackrefit + muonTable + muonExternalVecVarsTable + fsrTable)

# remove track refit stuff for low pu
run2_nanoAOD_LowPU.toReplaceWith(muonTables, muonTables.copyAndExclude([geopro, tracksfrommuons, trackrefit, muonExternalVecVarsTable]))
run2_nanoAOD_LowPU.toReplaceWith(muonMC, muonMC.copyAndExclude([trackrefitideal, mergedGlobalIdxs, muonIdealTable, muonIdealExternalVecVarsTable]))
run2_nanoAOD_LowPU.toModify(muonTable, externalVariables = cms.PSet())

run2_nanoAOD_LowPU.toReplaceWith(muonTables, muonTables.copyAndExclude([geopro, tracksfrommuons, trackrefit, trackrefitbs, mergedGlobalIdxs, muonExternalVecVarsTable]))

17 changes: 17 additions & 0 deletions PhysicsTools/NanoAOD/python/nano_cff.py
Original file line number Diff line number Diff line change
Expand Up @@ -443,13 +443,30 @@ def nanoAOD_customizeData(process):
process = nanoAOD_recalibrateMETs(process,isData=True)
for modifier in run2_nanoAOD_94XMiniAODv1, run2_nanoAOD_94XMiniAODv2:
modifier.toModify(process, lambda p: nanoAOD_runMETfixEE2017(p,isData=True))

# load geometry needed by g4e propagator
process.GlobalTag.toGet = cms.VPSet(
cms.PSet(
record = cms.string("GeometryFileRcd"),tag = cms.string("XMLFILE_Geometry_2016_81YV1_Extended2016_mc"),label = cms.untracked.string("Extended"),
),
)

# load 3d field map and use it for g4e propagator
from MagneticField.ParametrizedEngine.parametrizedMagneticField_PolyFit3D_cfi import ParametrizedMagneticFieldProducer as PolyFit3DMagneticFieldProducer
process.PolyFit3DMagneticFieldProducer = PolyFit3DMagneticFieldProducer
fieldlabel = "PolyFit3DMf"
process.PolyFit3DMagneticFieldProducer.label = fieldlabel
process.Geant4ePropagator.MagneticFieldLabel = fieldlabel

return process

def nanoAOD_customizeMC(process):
process = nanoAOD_customizeCommon(process)
process = nanoAOD_recalibrateMETs(process,isData=False)
for modifier in run2_nanoAOD_94XMiniAODv1, run2_nanoAOD_94XMiniAODv2:
modifier.toModify(process, lambda p: nanoAOD_runMETfixEE2017(p,isData=False))
# remove duplicate branch with incomplete content *TODO* make this more elegant
del process.muonExternalVecVarsTable.variables.cvhmergedGlobalIdxs
return process

###Customizations needed for Wmass analysis
Expand Down
8 changes: 6 additions & 2 deletions SimG4Core/GeometryProducer/interface/GeometryProducer.h
Original file line number Diff line number Diff line change
Expand Up @@ -26,6 +26,7 @@ class SimProducer;
class DDDWorld;
class G4RunManagerKernel;
class SimTrackManager;
class DDCompactView;

class GeometryProducer : public edm::one::EDProducer<edm::one::SharedResources, edm::one::WatchRuns> {
public:
Expand All @@ -46,9 +47,7 @@ class GeometryProducer : public edm::one::EDProducer<edm::one::SharedResources,
void updateMagneticField(edm::EventSetup const &es);

G4RunManagerKernel *m_kernel;
bool m_pUseMagneticField;
edm::ParameterSet m_pField;
bool m_pUseSensitiveDetectors;
SimActivityRegistry m_registry;
std::vector<std::shared_ptr<SimWatcher>> m_watchers;
std::vector<std::shared_ptr<SimProducer>> m_producers;
Expand All @@ -58,7 +57,12 @@ class GeometryProducer : public edm::one::EDProducer<edm::one::SharedResources,
std::vector<SensitiveTkDetector *> m_sensTkDets;
std::vector<SensitiveCaloDetector *> m_sensCaloDets;
edm::ParameterSet m_p;

mutable const DDCompactView *m_pDD;

bool m_firstRun;
bool m_pUseMagneticField;
bool m_pUseSensitiveDetectors;
};

#endif
Loading

0 comments on commit 9a6c375

Please sign in to comment.