Skip to content
Merged
Show file tree
Hide file tree
Changes from 1 commit
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -42,6 +42,7 @@ struct LtrCalibData {
std::vector<uint16_t> matchedLtrIDs; ///< matched laser track IDs
std::vector<uint16_t> nTrackTF; ///< number of laser tracks per TF
std::vector<float> dEdx; ///< dE/dx of each track
float tp{0.f}; ///< temperature over pressure ratio

bool isValid() const
{
Expand Down Expand Up @@ -138,7 +139,7 @@ struct LtrCalibData {
dEdx.clear();
}

ClassDefNV(LtrCalibData, 4);
ClassDefNV(LtrCalibData, 5);
};

} // namespace o2::tpc
Expand Down
29 changes: 19 additions & 10 deletions DataFormats/Detectors/TPC/include/DataFormatsTPC/VDriftCorrFact.h
Original file line number Diff line number Diff line change
Expand Up @@ -26,27 +26,36 @@ namespace o2::tpc
{

struct VDriftCorrFact {
long firstTime{}; ///< first time stamp of processed TFs
long lastTime{}; ///< last time stamp of processed TFs
long creationTime{}; ///< time of creation
float corrFact{1.0}; ///< drift velocity correction factor (multiplicative)
float corrFactErr{0.0}; ///< stat error of correction factor
float refVDrift{0.}; ///< reference vdrift for which factor was extracted
long firstTime{}; ///< first time stamp of processed TFs
long lastTime{}; ///< last time stamp of processed TFs
long creationTime{}; ///< time of creation
float corrFact{1.0}; ///< drift velocity correction factor (multiplicative)
float corrFactErr{0.0}; ///< stat error of correction factor
float refVDrift{0.}; ///< reference vdrift for which factor was extracted
float refTimeOffset{0.}; ///< additive time offset reference (\mus)
float timeOffsetCorr{0.}; ///< additive time offset correction (\mus)
float refTP{0.}; ///< reference temperature / pressure for which refVDrift was extracted

float getVDrift() const { return refVDrift * corrFact; }
float getVDriftError() const { return refVDrift * corrFactErr; }

float getTimeOffset() const { return refTimeOffset + timeOffsetCorr; }

// renormalize VDrift reference and correction either to provided new reference (if >0) or to correction 1 wrt current reference
void normalize(float newVRef = 0.f)
void normalize(float newVRef = 0.f, float tp = 0.f)
{
float normVDrift = newVRef;
if (newVRef == 0.f) {
newVRef = refVDrift * corrFact;
normVDrift = refVDrift * corrFact;
if ((tp == 0) || (refTP == 0)) {
newVRef = normVDrift; // no T/P scaling applied
} else {
// linear scaling based on relative change of T/P
newVRef = normVDrift * (1 + (tp - refTP) / refTP);
Copy link
Copy Markdown
Collaborator

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

@matthias-kleiner I am wondering if this is what we want. The original reason behind this complicated normalizations was to have
a reasonable refVDrift constant over the whole, so that the residuals extraction calibration, which have granularity different from that of VDrift, would use the same refVDrift value, see

// Attention! For the refit we are using reference VDrift and TDriftOffest rather than high-rate calibrated, since we want to have fixed reference over the run
if (v.refVDrift != mTPCVDriftRef) {
mTPCVDriftRef = v.refVDrift;
mTPCDriftTimeOffsetRef = v.refTimeOffset;
LOGP(info, "Imposing reference VDrift={}/TDrift={} for TPC residuals extraction", mTPCVDriftRef, mTPCDriftTimeOffsetRef);
o2::tpc::TPCFastTransformHelperO2::instance()->updateCalibration(*mFastTransform, 0, 1.0, mTPCVDriftRef, mTPCDriftTimeOffsetRef);
}
.
Therefore, the variation of the full VDrift was provided via time-dependent corrFact.

In the current code, the normalize method does not modify the result of getVDrift(), it is just

  1. setting the refVDrift to the requested newVRef if it is non-0 and readjusts the corrFact to have the full VDrift unchanged (used to make sure that the next calibration slot creates VDrift with the same reference),
    or
  2. at the 1st call w/o argument normalizes refVDrift to full VDrift and therefore sets the corrections to 1 (to have the refVDrift as correct as possible).

This change will not only scale the final getVDrift() value by the relative T/P change (which is OK) but will also modify the refVDrift according to the current meteo condition. Should not instead just corrFact (and refTP) adjusted to account for changed T/P?

Copy link
Copy Markdown
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

Hi @shahor02 , thanks for the comment and information. I was not aware of that specific reason why we wanted to have a constant refVDrift. In that case, we should indeed simply adjust the corrFact to account for changes in the T/P.

refTP = tp; // update reference T/P
}
}
float fact = refVDrift / newVRef;
float fact = refVDrift / normVDrift;
refVDrift = newVRef;
corrFactErr *= fact;
corrFact *= fact;
Expand All @@ -66,7 +75,7 @@ struct VDriftCorrFact {
}
}

ClassDefNV(VDriftCorrFact, 2);
ClassDefNV(VDriftCorrFact, 3);
};

} // namespace o2::tpc
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -331,6 +331,9 @@ struct TimeSeriesdEdx {
};

struct TimeSeriesITSTPC {
float mVDrift = 0; ///< drift velocity in cm/us
float mPressure = 0; ///< pressure
float mTemperature = 0; ///< temperature
TimeSeries mTSTPC; ///< TPC standalone DCAs
TimeSeries mTSITSTPC; ///< ITS-TPC standalone DCAs
ITSTPC_Matching mITSTPCAll; ///< ITS-TPC matching efficiency for ITS standalone + afterburner
Expand Down Expand Up @@ -499,7 +502,7 @@ struct TimeSeriesITSTPC {
nVertexContributors_Quantiles.resize(nTotalQ);
}

ClassDefNV(TimeSeriesITSTPC, 5);
ClassDefNV(TimeSeriesITSTPC, 6);
};

} // end namespace tpc
Expand Down
2 changes: 1 addition & 1 deletion Detectors/GlobalTracking/src/MatchTPCITS.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -101,7 +101,7 @@ void MatchTPCITS::run(const o2::globaltracking::RecoContainer& inp,
break;
}
if (mVDriftCalibOn) { // in the beginning of the output vector we send the full and reference VDrift used for this TF
calib.emplace_back(mTPCVDrift, mTPCDrift.refVDrift, -999.);
calib.emplace_back(mTPCVDrift, mTPCDrift.refVDrift, mTPCDrift.refTP);
calib.emplace_back(mTPCDriftTimeOffset, mTPCDrift.refTimeOffset, -999.);
}

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -25,7 +25,6 @@
#include <gsl/span>
#include <string_view>

#include "CommonConstants/MathConstants.h"
#include "CommonUtils/TreeStreamRedirector.h"
#include "DataFormatsTPC/TrackTPC.h"
#include "DataFormatsTPC/LaserTrack.h"
Expand Down Expand Up @@ -74,10 +73,12 @@ class CalibLaserTracks
~CalibLaserTracks() = default;

/// process all tracks of one TF
void fill(const gsl::span<const TrackTPC> tracks);
/// \param tp ratio of temperature over pressure
void fill(const gsl::span<const TrackTPC> tracks, float tp = 0);

/// process all tracks of one TF
void fill(std::vector<TrackTPC> const& tracks);
/// \param tp ratio of temperature over pressure
void fill(std::vector<TrackTPC> const& tracks, float tp = 0);

/// process single track
void processTrack(const TrackTPC& track);
Expand Down Expand Up @@ -163,6 +164,8 @@ class CalibLaserTracks
float mDriftV{0}; ///< drift velocity used during reconstruction
float mTOffsetMUS{0}; ///< time offset in \mus to impose
float mZbinWidth{0}; ///< width of a bin in us
float mAvgTP{0}; ///< ratio of average temperature over pressure
float mAvgDriftV{0}; ///< average drift velocity used for the laser track calibration
uint64_t mTFstart{0}; ///< start time of processed time frames
uint64_t mTFend{0}; ///< end time of processed time frames
LtrCalibData mCalibDataTF{}; ///< calibration data for single TF (debugging)
Expand All @@ -184,7 +187,7 @@ class CalibLaserTracks
/// perform fits on the matched z-position pairs to extract the drift velocity correction factor and trigger offset
void fillCalibData(LtrCalibData& calibData, const std::vector<TimePair>& pairsA, const std::vector<TimePair>& pairsC);

ClassDefNV(CalibLaserTracks, 1);
ClassDefNV(CalibLaserTracks, 2);
};

} // namespace o2::tpc
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -63,15 +63,35 @@ class PressureTemperatureHelper
/// get pressure for given time stamp in ms
float getPressure(const ULong64_t timestamp) const { return interpolate(mPressure.second, mPressure.first, timestamp); }

/// manually set the pressure
void setPressure(const std::pair<std::vector<float>, std::vector<ULong64_t>>& pressure) { mPressure = pressure; }

/// manually set the temperature
void setTemperature(const std::pair<std::vector<float>, std::vector<ULong64_t>>& temperatureA, const std::pair<std::vector<float>, std::vector<ULong64_t>>& temperatureC)
{
mTemperatureA = temperatureA;
mTemperatureC = temperatureC;
}

/// get temperature for given time stamp in ms
dataformats::Pair<float, float> getTemperature(const ULong64_t timestamp) const { return dataformats::Pair<float, float>{interpolate(mTemperatureA.second, mTemperatureA.first, timestamp), interpolate(mTemperatureC.second, mTemperatureC.first, timestamp)}; }

/// get mean temperature over A and C side
float getMeanTemperature(const ULong64_t timestamp) const;

// get ratio of temperature over pressure for given time stamp
float getTP(int64_t ts) const;

static constexpr o2::header::DataDescription getDataDescriptionPressure() { return o2::header::DataDescription{"pressure"}; }
static constexpr o2::header::DataDescription getDataDescriptionTemperature() { return o2::header::DataDescription{"temperature"}; }

/// get minimum and maximum time stamps of the pressure and temperature data
std::pair<ULong64_t, ULong64_t> getMinMaxTime() const;

protected:
static void addInput(std::vector<o2::framework::InputSpec>& inputs, o2::framework::InputSpec&& isp);
static void addOutput(std::vector<o2::framework::OutputSpec>& outputs, o2::framework::OutputSpec&& osp);
static constexpr float toKelvin(float celsius) { return celsius + 273.15f; } // convert Celsius to Kelvin

std::pair<std::vector<float>, std::vector<ULong64_t>> mPressure; ///< pressure values for both measurements
std::pair<std::vector<float>, std::vector<ULong64_t>> mTemperatureA; ///< temperature values A-side
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -32,6 +32,7 @@ struct TPCVDTglContainer {
double driftVFullMean = 0.;
static float tOffsetRef;
static float driftVRef;
float tp = 0;

TPCVDTglContainer(int ntgl, float tglMax, int ndtgl, float dtglMax)
{
Expand All @@ -42,9 +43,12 @@ struct TPCVDTglContainer {
{
histo = std::make_unique<o2::dataformats::FlatHisto2D_f>(*(src.histo.get()));
entries = src.entries;
tp = src.tp;
driftVFullMean = src.driftVFullMean;
}

void fill(const gsl::span<const o2::dataformats::Triplet<float, float, float>> data)
/// \param tp ratio of temperature over pressure
void fill(const gsl::span<const o2::dataformats::Triplet<float, float, float>> data, float currentTemperaturePressure = 0)
{
if (data.size() < 3) { // first 2 entres always contains the {full and reference VDrift} and {full and reference DriftTimeOffset} used for the TF
return;
Expand All @@ -59,12 +63,14 @@ struct TPCVDTglContainer {
}
//
float vfull = data[0].first, vref = data[0].second;
const float temperaturePressure = (data[0].third == 0) ? currentTemperaturePressure : data[0].third;
if (driftVRef == 0.f) {
driftVRef = vref;
} else if (driftVRef != vref) {
LOGP(warn, "data with VDriftRef={} were received while initially was set to {}, keep old one", vref, driftVRef);
}
driftVFullMean = (driftVFullMean * nTFProc + vfull) / (nTFProc + 1);
tp = (tp * nTFProc + temperaturePressure) / (nTFProc + 1);
if (tOffsetRef == 0.f) {
tOffsetRef = data[1].first; // assign 1st full toffset as a reference
}
Expand All @@ -73,6 +79,11 @@ struct TPCVDTglContainer {

void merge(const TPCVDTglContainer* other)
{
const int norm = nTFProc + other->nTFProc;
if (norm > 0) {
tp = (tp * nTFProc + other->tp * other->nTFProc) / norm;
driftVFullMean = (driftVFullMean * nTFProc + other->driftVFullMean * other->nTFProc) / norm;
}
entries += other->entries;
histo->add(*(other->histo));
LOGP(debug, "Old entries:{} New entries:{} oldSum: {} newSum: {}", other->entries, entries, other->histo->getSum(), histo->getSum());
Expand All @@ -82,7 +93,7 @@ struct TPCVDTglContainer {
{
LOG(info) << "Nentries = " << entries;
}
ClassDefNV(TPCVDTglContainer, 1);
ClassDefNV(TPCVDTglContainer, 2);
};

class TPCVDriftTglCalibration : public o2::calibration::TimeSlotCalibration<TPCVDTglContainer>
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@

#include "GPUCommonRtypes.h"
#include "DataFormatsTPC/VDriftCorrFact.h"
#include "TPCCalibration/PressureTemperatureHelper.h"
#include <array>
#include <vector>
#include <string_view>
Expand Down Expand Up @@ -56,23 +57,28 @@ class VDriftHelper
Source getSource() const { return mSource; }
static std::string_view getSourceName(Source s) { return SourceNames[s]; }
std::string_view getSourceName() const { return SourceNames[mSource]; }
const auto& getPTHelper() const { return mPTHelper; }

bool accountCCDBInputs(const o2::framework::ConcreteDataMatcher& matcher, void* obj);
void extractCCDBInputs(o2::framework::ProcessingContext& pc, bool laser = true, bool itstpcTgl = true);
static void requestCCDBInputs(std::vector<o2::framework::InputSpec>& inputs, bool laser = true, bool itstpcTgl = true);

protected:
static void addInput(std::vector<o2::framework::InputSpec>& inputs, o2::framework::InputSpec&& isp);
bool extractTPForVDrift(VDriftCorrFact& vdrift, int64_t tsStepMS = 100 * 1000);
VDriftCorrFact mVDLaser{};
VDriftCorrFact mVDTPCITSTgl{};
VDriftCorrFact mVD{};
Source mSource{Source::Param}; // update source
bool mUpdated = false; // signal update, must be reset once new value is fetched
Source mSource{Source::Param}; // update source
bool mUpdated = false; // signal update, must be reset once new value is fetched
bool mIsTPScalingPossible = false; // if T/P scaling is possible always perform the updating
bool mForceParamDrift = false; // enforce vdrift from gasParam
bool mForceParamOffset = false; // enforce offset from DetectorParam
bool mForceTPScaling = false; // enforce T/P scaling from gasParam (scaling disabled by negative T or P)
uint32_t mMayRenormSrc = 0xffffffff; // if starting VDrift correction != 1, we will renorm reference in such a way that initial correction is 1.0, flag per source
PressureTemperatureHelper mPTHelper; // helper to extract pressure and temperature from CCDB

ClassDefNV(VDriftHelper, 1);
ClassDefNV(VDriftHelper, 2);
};
} // namespace o2::tpc
#endif
27 changes: 23 additions & 4 deletions Detectors/TPC/calibration/src/CalibLaserTracks.cxx
Original file line number Diff line number Diff line change
Expand Up @@ -24,13 +24,13 @@
#include <string_view>

using namespace o2::tpc;
void CalibLaserTracks::fill(std::vector<TrackTPC> const& tracks)
void CalibLaserTracks::fill(std::vector<TrackTPC> const& tracks, float tp)
{
fill(gsl::span(tracks.data(), tracks.size()));
fill(gsl::span(tracks.data(), tracks.size()), tp);
}

//______________________________________________________________________________
void CalibLaserTracks::fill(const gsl::span<const TrackTPC> tracks)
void CalibLaserTracks::fill(const gsl::span<const TrackTPC> tracks, float tp)
{
// ===| clean up TF data |===
mZmatchPairsTFA.clear();
Expand Down Expand Up @@ -63,6 +63,9 @@ void CalibLaserTracks::fill(const gsl::span<const TrackTPC> tracks)
mCalibDataTF.firstTime = mTFstart;
mCalibDataTF.lastTime = tfEnd;

mAvgTP = (mAvgTP * mCalibData.processedTFs + tp) / (mCalibData.processedTFs + 1);
mAvgDriftV = (mAvgDriftV * mCalibData.processedTFs + mDriftV) / (mCalibData.processedTFs + 1);

// ===| TF counters |===
++mCalibData.processedTFs;
++mCalibDataTF.processedTFs;
Expand Down Expand Up @@ -147,6 +150,8 @@ void CalibLaserTracks::processTrack(const TrackTPC& track)
<< "ltr=" << ltr // matched ideal laser track
<< "trOutLtr=" << parOutAtLtr // track rotated and propagated to ideal track position
<< "TPCTracks=" << writeTrack // original TPC track
<< "mDriftV=" << mDriftV
<< "laserTrackID=" << laserTrackID
<< "\n";
}
}
Expand Down Expand Up @@ -277,6 +282,18 @@ void CalibLaserTracks::merge(const CalibLaserTracks* other)
mCalibData.firstTime = std::min(mCalibData.firstTime, other->mCalibData.firstTime);
mCalibData.lastTime = std::max(mCalibData.lastTime, other->mCalibData.lastTime);

if ((mAvgTP > 0) && (other->mAvgTP > 0)) {
mAvgTP = (mAvgTP + other->mAvgTP) / 2.0;
} else if (other->mAvgTP > 0) {
mAvgTP = other->mAvgTP;
}

if ((mAvgDriftV > 0) && (other->mAvgDriftV > 0)) {
mAvgDriftV = (mAvgDriftV + other->mAvgDriftV) / 2.0;
} else if (other->mAvgDriftV > 0) {
mAvgDriftV = other->mAvgDriftV;
}

sort(mZmatchPairsA);
sort(mZmatchPairsC);

Expand All @@ -296,6 +313,7 @@ void CalibLaserTracks::endTF()
<< "zPairsA=" << mZmatchPairsTFA
<< "zPairsC=" << mZmatchPairsTFC
<< "calibData=" << mCalibDataTF
<< "mDriftV=" << mDriftV
<< "\n";
}
}
Expand Down Expand Up @@ -330,7 +348,7 @@ void CalibLaserTracks::fillCalibData(LtrCalibData& calibData, const std::vector<
auto dvA = fit(pairsA, "A-Side");
auto dvC = fit(pairsC, "C-Side");
calibData.creationTime = std::chrono::duration_cast<std::chrono::milliseconds>(std::chrono::system_clock::now().time_since_epoch()).count();
calibData.refVDrift = mDriftV;
calibData.refVDrift = mAvgDriftV;
calibData.dvOffsetA = dvA.x1;
calibData.dvCorrectionA = dvA.x2;
calibData.nTracksA = uint16_t(pairsA.size());
Expand All @@ -340,6 +358,7 @@ void CalibLaserTracks::fillCalibData(LtrCalibData& calibData, const std::vector<
calibData.nTracksC = uint16_t(pairsC.size());

calibData.refTimeOffset = mTOffsetMUS;
calibData.tp = mAvgTP;
}

//______________________________________________________________________________
Expand Down
Loading