Skip to content

Commit

Permalink
Refactor Charge Reweighting code
Browse files Browse the repository at this point in the history
  • Loading branch information
@tvami
tvami committed Dec 16, 2022
1 parent 52c0dac commit 75bf0e8
Show file tree
Hide file tree
Showing 21 changed files with 133 additions and 1,163 deletions.
3 changes: 3 additions & 0 deletions SimTracker/Common/BuildFile.xml
View file
Expand Up @@ -4,6 +4,9 @@
<use name="SimDataFormats/TrackingAnalysis"/>
<use name="SimDataFormats/TrackingHit"/>
<use name="SimDataFormats/CrossingFrame"/>
<use name="CondFormats/SiPixelObjects"/>
<use name="CondFormats/SiPixelTransient"/>
<use name="gsl"/>
<export>
<lib name="1"/>
</export>
130 changes: 0 additions & 130 deletions SimTracker/SiPhase2Digitizer/plugins/DigitizerUtility.h
View file

This file was deleted.

2 changes: 1 addition & 1 deletion SimTracker/SiPhase2Digitizer/plugins/PSPDigitizerAlgorithm.cc
View file
Expand Up @@ -52,7 +52,7 @@ bool PSPDigitizerAlgorithm::select_hit(const PSimHit& hit, double tCorr, double&
//
// -- Compare Signal with Threshold
//
bool PSPDigitizerAlgorithm::isAboveThreshold(const DigitizerUtility::SimHitInfo* hitInfo,
bool PSPDigitizerAlgorithm::isAboveThreshold(const digitizerUtility::SimHitInfo* hitInfo,
float charge,
float thr) const {
return (charge >= thr);
Expand Down
2 changes: 1 addition & 1 deletion SimTracker/SiPhase2Digitizer/plugins/PSPDigitizerAlgorithm.h
View file
Expand Up @@ -15,7 +15,7 @@ class PSPDigitizerAlgorithm : public Phase2TrackerDigitizerAlgorithm {
void init(const edm::EventSetup& es) override;

bool select_hit(const PSimHit& hit, double tCorr, double& sigScale) const override;
bool isAboveThreshold(const DigitizerUtility::SimHitInfo* hitInfo, float charge, float thr) const override;
bool isAboveThreshold(const digitizerUtility::SimHitInfo* hitInfo, float charge, float thr) const override;
bool isInBiasRailRegion(const PSimHit& hit) const;
void module_killing_DB(const Phase2TrackerGeomDetUnit* pixdet) override;

Expand Down
2 changes: 1 addition & 1 deletion SimTracker/SiPhase2Digitizer/plugins/PSSDigitizerAlgorithm.cc
View file
Expand Up @@ -57,7 +57,7 @@ bool PSSDigitizerAlgorithm::select_hit(const PSimHit& hit, double tCorr, double&
//
// -- Compare Signal with Threshold
//
bool PSSDigitizerAlgorithm::isAboveThreshold(const DigitizerUtility::SimHitInfo* const hisInfo,
bool PSSDigitizerAlgorithm::isAboveThreshold(const digitizerUtility::SimHitInfo* const hisInfo,
float charge,
float thr) const {
return (charge >= thr);
Expand Down
2 changes: 1 addition & 1 deletion SimTracker/SiPhase2Digitizer/plugins/PSSDigitizerAlgorithm.h
View file
Expand Up @@ -16,7 +16,7 @@ class PSSDigitizerAlgorithm : public Phase2TrackerDigitizerAlgorithm {
void init(const edm::EventSetup& es) override;

bool select_hit(const PSimHit& hit, double tCorr, double& sigScale) const override;
bool isAboveThreshold(const DigitizerUtility::SimHitInfo* hitInfo, float charge, float thr) const override;
bool isAboveThreshold(const digitizerUtility::SimHitInfo* hitInfo, float charge, float thr) const override;
void module_killing_DB(const Phase2TrackerGeomDetUnit* pixdet) override;

private:
Expand Down
14 changes: 7 additions & 7 deletions SimTracker/SiPhase2Digitizer/plugins/Phase2TrackerDigitizer.cc
View file
Expand Up @@ -24,7 +24,7 @@
#include "SimTracker/SiPhase2Digitizer/plugins/PixelDigitizerAlgorithm.h"
#include "SimTracker/SiPhase2Digitizer/plugins/Pixel3DDigitizerAlgorithm.h"
#include "SimTracker/SiPhase2Digitizer/plugins/PixelBrickedDigitizerAlgorithm.h"
#include "SimTracker/SiPhase2Digitizer/plugins/DigitizerUtility.h"
#include "SimTracker/Common/interface/DigitizerUtility.h"

#include "FWCore/Framework/interface/EventSetup.h"
#include "FWCore/ParameterSet/interface/ParameterSet.h"
Expand Down Expand Up @@ -296,13 +296,13 @@ namespace cms {
algotype != AlgorithmType::InnerPixelBricked)
continue;

std::map<int, DigitizerUtility::DigiSimInfo> digi_map;
std::map<int, digitizerUtility::DigiSimInfo> digi_map;
fiter->second->digitize(dynamic_cast<const Phase2TrackerGeomDetUnit*>(det_u), digi_map, tTopo_);

edm::DetSet<PixelDigi> collector(rawId);
edm::DetSet<PixelDigiSimLink> linkcollector(rawId);
for (auto const& digi_p : digi_map) {
DigitizerUtility::DigiSimInfo info = digi_p.second;
digitizerUtility::DigiSimInfo info = digi_p.second;
const auto& ip = PixelDigi::channelToPixel(digi_p.first);
collector.data.emplace_back(ip.first, ip.second, info.sig_tot);
for (auto const& sim_p : info.simInfoList) {
Expand Down Expand Up @@ -331,7 +331,7 @@ namespace cms {
}
} // namespace cms
namespace {
void addToCollector(edm::DetSet<PixelDigi>& collector, const int channel, const DigitizerUtility::DigiSimInfo& info) {
void addToCollector(edm::DetSet<PixelDigi>& collector, const int channel, const digitizerUtility::DigiSimInfo& info) {
// For premixing stage1 the channel must be decoded with PixelDigi
// so that when the row and column are inserted to PixelDigi the
// coded channel stays the same (so that it can then be decoded
Expand All @@ -341,7 +341,7 @@ namespace {
}
void addToCollector(edm::DetSet<Phase2TrackerDigi>& collector,
const int channel,
const DigitizerUtility::DigiSimInfo& info) {
const digitizerUtility::DigiSimInfo& info) {
const auto& ip = Phase2TrackerDigi::channelToPixel(channel);
collector.data.emplace_back(ip.first, ip.second, info.ot_bit);
}
Expand All @@ -360,13 +360,13 @@ namespace cms {
algotype == AlgorithmType::InnerPixelBricked)
continue;

std::map<int, DigitizerUtility::DigiSimInfo> digi_map;
std::map<int, digitizerUtility::DigiSimInfo> digi_map;
fiter->second->digitize(dynamic_cast<const Phase2TrackerGeomDetUnit*>(det_u), digi_map, tTopo_);

edm::DetSet<DigiType> collector(rawId);
edm::DetSet<PixelDigiSimLink> linkcollector(rawId);
for (auto const& digi_p : digi_map) {
DigitizerUtility::DigiSimInfo info = digi_p.second;
digitizerUtility::DigiSimInfo info = digi_p.second;
addToCollector(collector, digi_p.first, info);
for (auto const& sim_p : info.simInfoList) {
linkcollector.data.emplace_back(digi_p.first,
Expand Down
40 changes: 20 additions & 20 deletions SimTracker/SiPhase2Digitizer/plugins/Phase2TrackerDigitizerAlgorithm.cc
View file
Expand Up @@ -213,7 +213,7 @@ void Phase2TrackerDigitizerAlgorithm::accumulateSimHits(std::vector<PSimHit>::co
// Divide the track into small sub-segments
//
// =================================================================
std::vector<DigitizerUtility::EnergyDepositUnit> Phase2TrackerDigitizerAlgorithm::primary_ionization(
std::vector<digitizerUtility::EnergyDepositUnit> Phase2TrackerDigitizerAlgorithm::primary_ionization(
const PSimHit& hit) const {
// Straight line approximation for trajectory inside active media
constexpr float SegmentLength = 0.0010; // in cm (10 microns)
Expand Down Expand Up @@ -241,15 +241,15 @@ std::vector<DigitizerUtility::EnergyDepositUnit> Phase2TrackerDigitizerAlgorithm
elossVector.resize(NumberOfSegments, averageEloss);
}

std::vector<DigitizerUtility::EnergyDepositUnit> ionization_points;
std::vector<digitizerUtility::EnergyDepositUnit> ionization_points;
ionization_points.reserve(NumberOfSegments); // set size
// loop over segments
for (size_t i = 0; i < elossVector.size(); ++i) {
// Divide the segment into equal length subsegments
Local3DPoint point = hit.entryPoint() + ((i + 0.5) / NumberOfSegments) * direction;
float energy = elossVector[i] / GeVperElectron_; // Convert charge to elec.

DigitizerUtility::EnergyDepositUnit edu(energy, point); // define position,energy point
digitizerUtility::EnergyDepositUnit edu(energy, point); // define position,energy point
ionization_points.push_back(edu); // save
LogDebug("Phase2TrackerDigitizerAlgorithm")
<< i << " " << edu.x() << " " << edu.y() << " " << edu.z() << " " << edu.energy();
Expand Down Expand Up @@ -319,14 +319,14 @@ std::vector<float> Phase2TrackerDigitizerAlgorithm::fluctuateEloss(
// Include the effect of E-field and B-field
//
// =====================================================================
std::vector<DigitizerUtility::SignalPoint> Phase2TrackerDigitizerAlgorithm::drift(
std::vector<digitizerUtility::SignalPoint> Phase2TrackerDigitizerAlgorithm::drift(
const PSimHit& hit,
const Phase2TrackerGeomDetUnit* pixdet,
const GlobalVector& bfield,
const std::vector<DigitizerUtility::EnergyDepositUnit>& ionization_points) const {
const std::vector<digitizerUtility::EnergyDepositUnit>& ionization_points) const {
LogDebug("Phase2TrackerDigitizerAlgorithm") << "enter drift ";

std::vector<DigitizerUtility::SignalPoint> collection_points;
std::vector<digitizerUtility::SignalPoint> collection_points;
collection_points.reserve(ionization_points.size()); // set size
LocalVector driftDir = DriftDirection(pixdet, bfield, hit.detUnitId()); // get the charge drift direction
if (driftDir.z() == 0.) {
Expand Down Expand Up @@ -406,7 +406,7 @@ std::vector<DigitizerUtility::SignalPoint> Phase2TrackerDigitizerAlgorithm::drif
LogDebug("Phase2TrackerDigitizerAlgorithm")
<< "Dift DistanceZ = " << driftDistance << " module thickness = " << moduleThickness
<< " Start Energy = " << val.energy() << " Energy after loss= " << energyOnCollector;
DigitizerUtility::SignalPoint sp(CloudCenterX, CloudCenterY, Sigma_x, Sigma_y, hit.tof(), energyOnCollector);
digitizerUtility::SignalPoint sp(CloudCenterX, CloudCenterY, Sigma_x, Sigma_y, hit.tof(), energyOnCollector);

// Load the Charge distribution parameters
collection_points.push_back(sp);
Expand All @@ -422,7 +422,7 @@ void Phase2TrackerDigitizerAlgorithm::induce_signal(
const size_t hitIndex,
const uint32_t tofBin,
const Phase2TrackerGeomDetUnit* pixdet,
const std::vector<DigitizerUtility::SignalPoint>& collection_points) {
const std::vector<digitizerUtility::SignalPoint>& collection_points) {
// X - Rows, Left-Right, 160, (1.6cm) for barrel
// Y - Columns, Down-Up, 416, (6.4cm)
const Phase2TrackerTopology* topol = &pixdet->specificTopology();
Expand Down Expand Up @@ -566,9 +566,9 @@ void Phase2TrackerDigitizerAlgorithm::induce_signal(
float corr_time = hit.tof() - pixdet->surface().toGlobal(hit.localPosition()).mag() * c_inv;
for (auto const& hit_s : hit_signal) {
int chan = hit_s.first;
theSignal[chan] +=
(makeDigiSimLinks_ ? DigitizerUtility::Amplitude(hit_s.second, &hit, hit_s.second, corr_time, hitIndex, tofBin)
: DigitizerUtility::Amplitude(hit_s.second, nullptr, hit_s.second));
theSignal[chan] += (makeDigiSimLinks_ ? digitizerUtility::Ph2Amplitude(
hit_s.second, &hit, hit_s.second, corr_time, hitIndex, tofBin)
: digitizerUtility::Ph2Amplitude(hit_s.second, nullptr, hit_s.second));
}
}
// ======================================================================
Expand Down Expand Up @@ -616,13 +616,13 @@ void Phase2TrackerDigitizerAlgorithm::add_cross_talk(const Phase2TrackerGeomDetU
auto XtalkPrev = std::make_pair(hitChan.first - 1, hitChan.second);
int chanXtalkPrev = pixelFlag_ ? PixelDigi::pixelToChannel(XtalkPrev.first, XtalkPrev.second)
: Phase2TrackerDigi::pixelToChannel(XtalkPrev.first, XtalkPrev.second);
signalNew.emplace(chanXtalkPrev, DigitizerUtility::Amplitude(signalInElectrons_Xtalk, nullptr, -1.0));
signalNew.emplace(chanXtalkPrev, digitizerUtility::Ph2Amplitude(signalInElectrons_Xtalk, nullptr, -1.0));
}
if (hitChan.first < numRows - 1) {
auto XtalkNext = std::make_pair(hitChan.first + 1, hitChan.second);
int chanXtalkNext = pixelFlag_ ? PixelDigi::pixelToChannel(XtalkNext.first, XtalkNext.second)
: Phase2TrackerDigi::pixelToChannel(XtalkNext.first, XtalkNext.second);
signalNew.emplace(chanXtalkNext, DigitizerUtility::Amplitude(signalInElectrons_Xtalk, nullptr, -1.0));
signalNew.emplace(chanXtalkNext, digitizerUtility::Ph2Amplitude(signalInElectrons_Xtalk, nullptr, -1.0));
}
}
for (auto const& l : signalNew) {
Expand All @@ -631,7 +631,7 @@ void Phase2TrackerDigitizerAlgorithm::add_cross_talk(const Phase2TrackerGeomDetU
if (iter != theSignal.end()) {
theSignal[chan] += l.second.ampl();
} else {
theSignal.emplace(chan, DigitizerUtility::Amplitude(l.second.ampl(), nullptr, -1.0));
theSignal.emplace(chan, digitizerUtility::Ph2Amplitude(l.second.ampl(), nullptr, -1.0));
}
}
}
Expand Down Expand Up @@ -679,7 +679,7 @@ void Phase2TrackerDigitizerAlgorithm::add_noisy_cells(const Phase2TrackerGeomDet
<< " Storing noise = " << el.first << " " << el.second << " " << ix << " " << iy << " " << chan;

if (theSignal[chan] == 0)
theSignal[chan] = DigitizerUtility::Amplitude(el.second, nullptr, -1.);
theSignal[chan] = digitizerUtility::Ph2Amplitude(el.second, nullptr, -1.);
}
}
// ============================================================================
Expand Down Expand Up @@ -854,15 +854,15 @@ void Phase2TrackerDigitizerAlgorithm::loadAccumulator(uint32_t detId, const std:
// the input channel is always with PixelDigi definition
// if needed, that has to be converted to Phase2TrackerDigi convention
for (const auto& elem : accumulator) {
auto inserted = theSignal.emplace(elem.first, DigitizerUtility::Amplitude(elem.second, nullptr));
auto inserted = theSignal.emplace(elem.first, digitizerUtility::Ph2Amplitude(elem.second, nullptr));
if (!inserted.second) {
throw cms::Exception("LogicError") << "Signal was already set for DetId " << detId;
}
}
}

void Phase2TrackerDigitizerAlgorithm::digitize(const Phase2TrackerGeomDetUnit* pixdet,
std::map<int, DigitizerUtility::DigiSimInfo>& digi_map,
std::map<int, digitizerUtility::DigiSimInfo>& digi_map,
const TrackerTopology* tTopo) {
uint32_t detID = pixdet->geographicalId().rawId();
auto it = _signal.find(detID);
Expand Down Expand Up @@ -914,16 +914,16 @@ void Phase2TrackerDigitizerAlgorithm::digitize(const Phase2TrackerGeomDetUnit* p

// Digitize if the signal is greater than threshold
for (auto const& s : theSignal) {
const DigitizerUtility::Amplitude& sig_data = s.second;
const digitizerUtility::Ph2Amplitude& sig_data = s.second;
float signalInElectrons = sig_data.ampl();

const auto& info_list = sig_data.simInfoList();
const DigitizerUtility::SimHitInfo* hitInfo = nullptr;
const digitizerUtility::SimHitInfo* hitInfo = nullptr;
if (!info_list.empty())
hitInfo = std::max_element(info_list.begin(), info_list.end())->second.get();

if (isAboveThreshold(hitInfo, signalInElectrons, theThresholdInE)) { // check threshold
DigitizerUtility::DigiSimInfo info;
digitizerUtility::DigiSimInfo info;
info.sig_tot = convertSignalToAdc(detID, signalInElectrons, theThresholdInE); // adc
info.ot_bit = signalInElectrons > theHIPThresholdInE ? true : false;
if (makeDigiSimLinks_) {
Expand Down

0 comments on commit 75bf0e8

Please sign in to comment.