From 4a2eae2fc593225d974cfd29ec0f55a1c3a5e18d Mon Sep 17 00:00:00 2001 From: shahoian Date: Fri, 12 Jul 2024 14:31:21 +0200 Subject: [PATCH 1/2] POD version of TPCFastSpaceChargeCorrection Created from the TPCFastTransform using: std::vector buff; // can be also pmr::vector from DPL make(..) const o2::gpu::TPCFastTransformPOD& pod = o2::gpu::TPCFastTransformPOD::create(buff, tpc_fast_transform); The buff vector will be expanded during creation and can be sent over DPL. On the receiving side, it should be cast as: const auto& podTransform = o2::gpu::TPCFastTransformPOD::get(pc.inputs().get>(ref)); No initialization is needed, the transform methods (at the moment all methods of TPCFastSpaceChargeCorrection are implemented) can be directly queried from the object received over sh.memory. Add perfmormance test for pod map --- GPU/TPCFastTransformation/CMakeLists.txt | 1 + .../TPCFastSpaceChargeCorrection.h | 2 + GPU/TPCFastTransformation/TPCFastTransform.h | 3 + .../TPCFastTransformPOD.cxx | 243 ++++++ .../TPCFastTransformPOD.h | 807 ++++++++++++++++++ .../TPCFastTransformationLinkDef_O2.h | 1 + 6 files changed, 1057 insertions(+) create mode 100644 GPU/TPCFastTransformation/TPCFastTransformPOD.cxx create mode 100644 GPU/TPCFastTransformation/TPCFastTransformPOD.h diff --git a/GPU/TPCFastTransformation/CMakeLists.txt b/GPU/TPCFastTransformation/CMakeLists.txt index b338e1492cc6c..c1c155805117e 100644 --- a/GPU/TPCFastTransformation/CMakeLists.txt +++ b/GPU/TPCFastTransformation/CMakeLists.txt @@ -25,6 +25,7 @@ set(SRCS TPCFastSpaceChargeCorrection.cxx TPCFastSpaceChargeCorrectionMap.cxx TPCFastTransform.cxx + TPCFastTransformPOD.cxx CorrectionMapsHelper.cxx ) diff --git a/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h b/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h index ecbc3c4ad73f5..5560b7d8d8fef 100644 --- a/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h +++ b/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h @@ -37,6 +37,8 @@ namespace gpu /// class TPCFastSpaceChargeCorrection : public FlatObject { + friend class TPCFastTransformPOD; + public: /// /// \brief The struct contains necessary info for TPC padrow diff --git a/GPU/TPCFastTransformation/TPCFastTransform.h b/GPU/TPCFastTransformation/TPCFastTransform.h index fe52bffe14acc..492702a185c08 100644 --- a/GPU/TPCFastTransformation/TPCFastTransform.h +++ b/GPU/TPCFastTransformation/TPCFastTransform.h @@ -250,6 +250,9 @@ class TPCFastTransform : public FlatObject /// Return TOF correction (vdrift / C) GPUd() float getTOFCorr() const { return mLdriftCorr; } + /// Return nominal PV Z position correction (vdrift / C) + GPUd() float getPrimVtxZ() const { return mPrimVtxZ; } + /// Return map lumi GPUd() float getLumi() const { return mLumi; } diff --git a/GPU/TPCFastTransformation/TPCFastTransformPOD.cxx b/GPU/TPCFastTransformation/TPCFastTransformPOD.cxx new file mode 100644 index 0000000000000..2435bcfb52e0d --- /dev/null +++ b/GPU/TPCFastTransformation/TPCFastTransformPOD.cxx @@ -0,0 +1,243 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file TPCFastTransformPOD.cxx +/// \brief Implementation of POD correction map +/// +/// \author ruben.shahoayn@cern.ch + +#include "TPCFastTransformPOD.h" +#include "GPUDebugStreamer.h" +#if !defined(GPUCA_GPUCODE) +#include +#endif + +namespace GPUCA_NAMESPACE +{ +namespace gpu +{ + +#if !defined(GPUCA_GPUCODE) + +size_t TPCFastTransformPOD::estimateSize(const TPCFastSpaceChargeCorrection& origCorr) +{ + // estimate size of own buffer + const size_t selfSizeFix = sizeof(TPCFastTransformPOD); + size_t nextDynOffs = alignOffset(selfSizeFix); + nextDynOffs = alignOffset(nextDynOffs + origCorr.mNumberOfScenarios * sizeof(size_t)); // spline scenarios start here + // space for splines + for (int isc = 0; isc < origCorr.mNumberOfScenarios; isc++) { + const auto& spline = origCorr.mScenarioPtr[isc]; + nextDynOffs = alignOffset(nextDynOffs + sizeof(spline)); + } + // space for splines data + for (int is = 0; is < 3; is++) { + for (int slice = 0; slice < origCorr.mGeo.getNumberOfSlices(); slice++) { + for (int row = 0; row < NROWS; row++) { + const auto& spline = origCorr.getSpline(slice, row); + int nPar = spline.getNumberOfParameters(); + if (is == 1) { + nPar = nPar / 3; + } + if (is == 2) { + nPar = nPar * 2 / 3; + } + nextDynOffs += nPar * sizeof(float); + } + } + } + nextDynOffs = alignOffset(nextDynOffs); + return nextDynOffs; +} + +TPCFastTransformPOD& TPCFastTransformPOD::create(char* buff, size_t buffSize, const TPCFastSpaceChargeCorrection& origCorr) +{ + // instantiate object to already created buffer of the right size + assert(buffSize > sizeof(TPCFastTransformPOD)); + auto& podMap = getNonConst(buff); + podMap.mApplyCorrections = true; // by default always apply corrections + + // copy fixed size data --- start + podMap.mNumberOfScenarios = origCorr.mNumberOfScenarios; + std::memcpy(&podMap.mGeo, &origCorr.mGeo, sizeof(TPCFastTransformGeo)); // copy geometry (fixed size) + for (int row = 0; row < NROWS; row++) { + podMap.mRowInfo[row] = origCorr.getRowInfo(row); // dataOffsetBytes will be modified later + } + for (int slice = 0; slice < TPCFastTransformGeo::getNumberOfSlices(); slice++) { + podMap.mSliceInfo[slice] = origCorr.getSliceInfo(slice); + for (int row = 0; row < NROWS; row++) { + podMap.mSliceRowInfo[NROWS * slice + row] = origCorr.getSliceRowInfo(slice, row); + } + } + podMap.mInterpolationSafetyMargin = origCorr.fInterpolationSafetyMargin; + podMap.mTimeStamp = origCorr.mTimeStamp; + // + // init data members coming from the TPCFastTrasform + podMap.mVdrift = 0.; + podMap.mT0 = 0.; + podMap.mVdriftCorrY = 0.; + podMap.mLdriftCorr = 0.; + podMap.mTOFcorr = 0.; + podMap.mPrimVtxZ = 0.; + // copy fixed size data --- end + + size_t nextDynOffs = alignOffset(sizeof(TPCFastTransformPOD)); + // copy slice scenarios + podMap.mOffsScenariosOffsets = nextDynOffs; // spline scenarios offsets start here + LOGP(debug, "Set mOffsScenariosOffsets = {}", podMap.mOffsScenariosOffsets); + nextDynOffs = alignOffset(nextDynOffs + podMap.mNumberOfScenarios * sizeof(size_t)); // spline scenarios start here + + // copy spline objects + size_t* scenOffs = reinterpret_cast(buff + podMap.mOffsScenariosOffsets); + for (int isc = 0; isc < origCorr.mNumberOfScenarios; isc++) { + scenOffs[isc] = nextDynOffs; + const auto& spline = origCorr.mScenarioPtr[isc]; + if (buffSize < nextDynOffs + sizeof(spline)) { + throw std::runtime_error(fmt::format("attempt to copy {} bytes for spline for scenario {} to {}, overflowing the buffer of size {}", sizeof(spline), isc, nextDynOffs + sizeof(spline), buffSize)); + } + std::memcpy(buff + scenOffs[isc], &spline, sizeof(spline)); + nextDynOffs = alignOffset(nextDynOffs + sizeof(spline)); + LOGP(debug, "Copy {} bytes for spline scenario {} (ptr:{}) to offsset {}", sizeof(spline), isc, (void*)&spline, scenOffs[isc]); + } + + // copy splines data + for (int is = 0; is < 3; is++) { + float* data = reinterpret_cast(buff + nextDynOffs); + LOGP(debug, "splinID={} start offset {} -> {}", is, nextDynOffs, (void*)data); + for (int slice = 0; slice < origCorr.mGeo.getNumberOfSlices(); slice++) { + podMap.mSplineDataOffsets[slice][is] = nextDynOffs; + size_t rowDataOffs = 0; + for (int row = 0; row < NROWS; row++) { + const auto& spline = origCorr.getSpline(slice, row); + const float* dataOr = origCorr.getSplineData(slice, row, is); + int nPar = spline.getNumberOfParameters(); + if (is == 1) { + nPar = nPar / 3; + } + if (is == 2) { + nPar = nPar * 2 / 3; + } + LOGP(debug, "Copying {} floats for spline{} of slice:{} row:{} to offset {}", nPar, is, slice, row, nextDynOffs); + size_t nbcopy = nPar * sizeof(float); + if (buffSize < nextDynOffs + nbcopy) { + throw std::runtime_error(fmt::format("attempt to copy {} bytes of data for spline{} of slice{}/row{} to {}, overflowing the buffer of size {}", nbcopy, is, slice, row, nextDynOffs, buffSize)); + } + std::memcpy(data, dataOr, nbcopy); + podMap.getRowInfo(row).dataOffsetBytes[is] = rowDataOffs; + rowDataOffs += nbcopy; + data += nPar; + nextDynOffs += nbcopy; + } + } + } + podMap.mTotalSize = alignOffset(nextDynOffs); + if (buffSize != podMap.mTotalSize) { + throw std::runtime_error(fmt::format("Estimated buffer size {} differs from filled one {}", buffSize, podMap.mTotalSize)); + } + return podMap; +} + +TPCFastTransformPOD& TPCFastTransformPOD::create(char* buff, size_t buffSize, const TPCFastTransform& src) +{ + // instantiate objec to already created buffer of the right size + auto& podMap = create(buff, buffSize, src.getCorrection()); + // set data members of TPCFastTransform + podMap.mVdrift = src.getVDrift(); + podMap.mT0 = src.getT0(); + podMap.mVdriftCorrY = src.getVdriftCorrY(); + podMap.mLdriftCorr = src.getLdriftCorr(); + podMap.mTOFcorr = src.getTOFCorr(); + podMap.mPrimVtxZ = src.getPrimVtxZ(); + // copy fixed size data --- end + return podMap; +} + +bool TPCFastTransformPOD::test(const TPCFastSpaceChargeCorrection& origCorr, int npoints) const +{ + if (npoints < 1) { + return false; + } + std::vector slice, row; + std::vector u, v, dxO, duO, dvO, dxP, duP, dvP, corrXO, corrXP, nomUO, nomVO, nomUP, nomVP; + slice.reserve(npoints); + row.reserve(npoints); + u.reserve(npoints); + v.reserve(npoints); + dxO.resize(npoints); + duO.resize(npoints); + dvO.resize(npoints); + corrXO.resize(npoints); + nomUO.resize(npoints); + nomVO.resize(npoints); + dxP.resize(npoints); + duP.resize(npoints); + dvP.resize(npoints); + corrXP.resize(npoints); + nomUP.resize(npoints); + nomVP.resize(npoints); + + for (int i = 0; i < npoints; i++) { + slice.push_back(gRandom->Integer(NSLICES)); + row.push_back(gRandom->Integer(NROWS)); + u.push_back(gRandom->Rndm() * 15); + v.push_back(gRandom->Rndm() * 200); + } + long origStart[3], origEnd[3], thisStart[3], thisEnd[3]; + origStart[0] = std::chrono::time_point_cast(std::chrono::system_clock::now()).time_since_epoch().count(); + for (int i = 0; i < npoints; i++) { + origCorr.getCorrection(slice[i], row[i], u[i], v[i], dxO[i], duO[i], dvO[i]); + } + origEnd[0] = origStart[1] = std::chrono::time_point_cast(std::chrono::system_clock::now()).time_since_epoch().count(); + for (int i = 0; i < npoints; i++) { + origCorr.getCorrectionInvCorrectedX(slice[i], row[i], u[i], v[i], corrXO[i]); + } + origEnd[1] = origStart[2] = std::chrono::time_point_cast(std::chrono::system_clock::now()).time_since_epoch().count(); + for (int i = 0; i < npoints; i++) { + origCorr.getCorrectionInvUV(slice[i], row[i], u[i], v[i], nomUO[i], nomVO[i]); + } + origEnd[2] = thisStart[0] = std::chrono::time_point_cast(std::chrono::system_clock::now()).time_since_epoch().count(); + for (int i = 0; i < npoints; i++) { + this->getCorrection(slice[i], row[i], u[i], v[i], dxP[i], duP[i], dvP[i]); + } + thisEnd[0] = thisStart[1] = std::chrono::time_point_cast(std::chrono::system_clock::now()).time_since_epoch().count(); + for (int i = 0; i < npoints; i++) { + this->getCorrectionInvCorrectedX(slice[i], row[i], u[i], v[i], corrXP[i]); + } + thisEnd[1] = thisStart[2] = std::chrono::time_point_cast(std::chrono::system_clock::now()).time_since_epoch().count(); + for (int i = 0; i < npoints; i++) { + this->getCorrectionInvUV(slice[i], row[i], u[i], v[i], nomUP[i], nomVP[i]); + } + thisEnd[2] = std::chrono::time_point_cast(std::chrono::system_clock::now()).time_since_epoch().count(); + // + size_t ndiff[3] = {}; + for (int i = 0; i < npoints; i++) { + if (dxO[i] != dxP[i] || duO[i] != duP[i] || dvO[i] != dvP[i]) { + ndiff[0]++; + } + if (corrXO[i] != corrXP[i]) { + ndiff[1]++; + } + if (nomUO[i] != nomUP[i] || nomVO[i] != nomVP[i]) { + ndiff[2]++; + } + } + // + LOGP(info, " (ns per call) original this Nmissmatch"); + LOGP(info, "getCorrection {:.3e} {:.3e} {}", double(origEnd[0] - origStart[0]) / npoints * 1000., double(thisEnd[0] - thisStart[0]) / npoints * 1000., ndiff[0]); + LOGP(info, "getCorrectionInvCorrectedX {:.3e} {:.3e} {}", double(origEnd[1] - origStart[1]) / npoints * 1000., double(thisEnd[1] - thisStart[1]) / npoints * 1000., ndiff[1]); + LOGP(info, "getCorrectionInvUV {:.3e} {:.3e} {}", double(origEnd[2] - origStart[2]) / npoints * 1000., double(thisEnd[2] - thisStart[2]) / npoints * 1000., ndiff[2]); + return ndiff[0] == 0 && ndiff[1] == 0 && ndiff[2] == 0; +} + +#endif + +} // namespace gpu +} // namespace GPUCA_NAMESPACE diff --git a/GPU/TPCFastTransformation/TPCFastTransformPOD.h b/GPU/TPCFastTransformation/TPCFastTransformPOD.h new file mode 100644 index 0000000000000..dd7e42b000c5a --- /dev/null +++ b/GPU/TPCFastTransformation/TPCFastTransformPOD.h @@ -0,0 +1,807 @@ +// Copyright 2019-2020 CERN and copyright holders of ALICE O2. +// See https://alice-o2.web.cern.ch/copyright for details of the copyright holders. +// All rights not expressly granted are reserved. +// +// This software is distributed under the terms of the GNU General Public +// License v3 (GPL Version 3), copied verbatim in the file "COPYING". +// +// In applying this license CERN does not waive the privileges and immunities +// granted to it by virtue of its status as an Intergovernmental Organization +// or submit itself to any jurisdiction. + +/// \file TPCFastTransformPOD.h +/// \brief POD correction map +/// +/// \author ruben.shahoayn@cern.ch + +#ifndef ALICEO2_GPU_TPCFastTransformPOD_H +#define ALICEO2_GPU_TPCFastTransformPOD_H + +#include "TPCFastTransform.h" + +/* +Binary buffer should be cast to TPCFastTransformPOD class using static TPCFastTransformPOD& t = get(buffer); method, +so that the its head becomes `this` pointer of the object. + +First we have all the fixed size data members mentioned explicitly. Part of them is duplicating fixed size +data members of TPCFastSpaceChargeCorrection but those starting with mOffs... provide the offset in bytes +(wrt this) for dynamic data which cannot be declared as data member explicitly (since we cannot have any +pointer except `this`) but obtained via getters using stored offsets wrt `this`. +This is followed dynamic part itself. + +dynamic part layout: +1) size_t[ mNumberOfScenarios ] array starting at offset mOffsScenariosOffsets, each element is the offset +of distict spline object (scenario in TPCFastSpaceChargeCorrection) +2) size_t[ mNSplineIDs ] array starting at offset mOffsSplineDataOffsets, each element is the offset of the +beginning of splines data for give splineID + +*/ + +namespace GPUCA_NAMESPACE +{ +namespace gpu +{ + +class TPCFastTransformPOD +{ + public: + using SplineType = TPCFastSpaceChargeCorrection::SplineType; + using RowInfo = TPCFastSpaceChargeCorrection::RowInfo; + using RowActiveArea = TPCFastSpaceChargeCorrection::RowActiveArea; + using SliceRowInfo = TPCFastSpaceChargeCorrection::SliceRowInfo; + using SliceInfo = TPCFastSpaceChargeCorrection::SliceInfo; + + /// convert prefilled buffer to TPCFastTransformPOD + GPUd() static const TPCFastTransformPOD& get(const char* head) { return *reinterpret_cast(head); } + + GPUd() int getApplyCorrections() const { return mApplyCorrections; } + GPUd() void setApplyCorrections(bool v = true) { mApplyCorrections = v; } + + /// _______________ high level methods a la TPCFastTransform _______________________ + /// + GPUd() void Transform(int slice, int row, float pad, float time, float& x, float& y, float& z, float vertexTime = 0) const; + GPUd() void TransformXYZ(int slice, int row, float& x, float& y, float& z) const; + + /// Transformation in the time frame + GPUd() void TransformInTimeFrame(int slice, int row, float pad, float time, float& x, float& y, float& z, float maxTimeBin) const; + GPUd() void TransformInTimeFrame(int slice, float time, float& z, float maxTimeBin) const; + + /// Inverse transformation + GPUd() void InverseTransformInTimeFrame(int slice, int row, float /*x*/, float y, float z, float& pad, float& time, float maxTimeBin) const; + + /// Inverse transformation: Transformed Y and Z -> transformed X + GPUd() void InverseTransformYZtoX(int slice, int row, float y, float z, float& x) const; + + /// Inverse transformation: Transformed Y and Z -> Y and Z, transformed w/o space charge correction + GPUd() void InverseTransformYZtoNominalYZ(int slice, int row, float y, float z, float& ny, float& nz) const; + + /// Inverse transformation: Transformed X, Y and Z -> X, Y and Z, transformed w/o space charge correction + GPUd() void InverseTransformXYZtoNominalXYZ(int slice, int row, float x, float y, float z, float& nx, float& ny, float& nz) const; + + /// Ideal transformation with Vdrift only - without calibration + GPUd() void TransformIdeal(int slice, int row, float pad, float time, float& x, float& y, float& z, float vertexTime) const; + GPUd() void TransformIdealZ(int slice, float time, float& z, float vertexTime) const; + + GPUd() void convPadTimeToUV(int slice, int row, float pad, float time, float& u, float& v, float vertexTime) const; + GPUd() void convPadTimeToUVinTimeFrame(int slice, int row, float pad, float time, float& u, float& v, float maxTimeBin) const; + GPUd() void convTimeToVinTimeFrame(int slice, float time, float& v, float maxTimeBin) const; + + GPUd() void convUVtoPadTime(int slice, int row, float u, float v, float& pad, float& time, float vertexTime) const; + GPUd() void convUVtoPadTimeInTimeFrame(int slice, int row, float u, float v, float& pad, float& time, float maxTimeBin) const; + GPUd() void convVtoTime(float v, float& time, float vertexTime) const; + + GPUd() float convTimeToZinTimeFrame(int slice, float time, float maxTimeBin) const; + GPUd() float convZtoTimeInTimeFrame(int slice, float z, float maxTimeBin) const; + GPUd() float convDeltaTimeToDeltaZinTimeFrame(int slice, float deltaTime) const; + GPUd() float convDeltaZtoDeltaTimeInTimeFrame(int slice, float deltaZ) const; + GPUd() float convDeltaZtoDeltaTimeInTimeFrameAbs(float deltaZ) const; + GPUd() float convZOffsetToVertexTime(int slice, float zOffset, float maxTimeBin) const; + GPUd() float convVertexTimeToZOffset(int slice, float vertexTime, float maxTimeBin) const; + + GPUd() void getTOFcorrection(int slice, int row, float x, float y, float z, float& dz) const; + + /// _______________ methods a la TPCFastSpaceChargeCorrection: cluster correction _______________________ + /// + GPUd() int getCorrection(int slice, int row, float u, float v, float& dx, float& du, float& dv) const; + + /// temporary method with the an way of calculating 2D spline + GPUd() int getCorrectionOld(int slice, int row, float u, float v, float& dx, float& du, float& dv) const; + + /// inverse correction: Corrected U and V -> coorrected X + GPUd() void getCorrectionInvCorrectedX(int slice, int row, float corrU, float corrV, float& corrX) const; + + /// inverse correction: Corrected U and V -> uncorrected U and V + GPUd() void getCorrectionInvUV(int slice, int row, float corrU, float corrV, float& nomU, float& nomV) const; + + /// maximal possible drift length of the active area + GPUd() float getMaxDriftLength(int slice, int row, float pad) const; + + /// maximal possible drift length of the active area + GPUd() float getMaxDriftLength(int slice, int row) const { return getSliceRowInfo(slice, row).activeArea.vMax; } + + /// maximal possible drift length of the active area + GPUd() float getMaxDriftLength(int slice) const { return getSliceInfo(slice).vMax; } + + /// maximal possible drift time of the active area + GPUd() float getMaxDriftTime(int slice, int row, float pad) const; + + /// maximal possible drift time of the active area + GPUd() float getMaxDriftTime(int slice, int row) const; + + /// maximal possible drift time of the active area + GPUd() float getMaxDriftTime(int slice) const; + + /// _______________ Utilities _______________________________________________ + + /// shrink u,v coordinats to the TPC row area +/- fkInterpolationSafetyMargin + GPUd() void schrinkUV(int slice, int row, float& u, float& v) const; + + /// shrink corrected u,v coordinats to the TPC row area +/- fkInterpolationSafetyMargin + GPUd() void schrinkCorrectedUV(int slice, int row, float& corrU, float& corrV) const; + + /// convert u,v to internal grid coordinates + GPUd() void convUVtoGrid(int slice, int row, float u, float v, float& gridU, float& gridV) const; + + /// convert u,v to internal grid coordinates + GPUd() void convGridToUV(int slice, int row, float gridU, float gridV, float& u, float& v) const; + + /// convert corrected u,v to internal grid coordinates + GPUd() void convCorrectedUVtoGrid(int slice, int row, float cu, float cv, float& gridU, float& gridV) const; + + /// TPC geometry information + GPUd() const TPCFastTransformGeo& getGeometry() const { return mGeo; } + + /// Gives its own size including dynamic part + GPUd() size_t size() const { return mTotalSize; } + + /// Gives the time stamp of the current calibaration parameters + GPUd() long int getTimeStamp() const { return mTimeStamp; } + + /// Gives the interpolation safety marging around the TPC row. + GPUd() float getInterpolationSafetyMargin() const { return mInterpolationSafetyMargin; } + + /// Return mVDrift in cm / time bin + GPUd() float getVDrift() const { return mVdrift; } + + /// Return T0 in time bin units + GPUd() float getT0() const { return mT0; } + + /// Return VdriftCorrY in time_bin / cn + GPUd() float getVdriftCorrY() const { return mVdriftCorrY; } + + /// Return LdriftCorr offset in cm + GPUd() float getLdriftCorr() const { return mLdriftCorr; } + + /// Return TOF correction + GPUd() float getTOFCorr() const { return mTOFcorr; } + + /// Return nominal PV Z position + GPUd() float getPrimVtxZ() const { return mPrimVtxZ; } + + /// Sets the time stamp of the current calibaration + GPUd() void setTimeStamp(long int v) { mTimeStamp = v; } + + /// Set safety marging for the interpolation around the TPC row. + /// Outside of this area the interpolation returns the boundary values. + GPUd() void setInterpolationSafetyMargin(float val) { mInterpolationSafetyMargin = val; } + + /// Gives TPC row info + GPUd() const RowInfo& getRowInfo(int row) const { return mRowInfo[row]; } + + /// Gives TPC slice info + GPUd() const SliceInfo& getSliceInfo(int slice) const { return mSliceInfo[slice]; } + + /// Gives a reference to a spline + GPUd() const SplineType& getSpline(int slice, int row) const { return *reinterpret_cast(getThis() + getScenarioOffset(getRowInfo(row).splineScenarioID)); } + + /// Gives pointer to spline data + GPUd() const float* getSplineData(int slice, int row, int iSpline = 0) const { return reinterpret_cast(getThis() + mSplineDataOffsets[slice][iSpline] + getRowInfo(row).dataOffsetBytes[iSpline]); } + + /// Gives TPC slice & row info + GPUd() const SliceRowInfo& getSliceRowInfo(int slice, int row) const { return mSliceRowInfo[NROWS * slice + row]; } + +#if !defined(GPUCA_GPUCODE) + /// Create POD transform from old flat-buffer one. Provided vector will serve as a buffer + template + static TPCFastTransformPOD& create(V& destVector, const TPCFastTransform& src); + + /// create filling only part corresponding to TPCFastSpaceChargeCorrection. Data members coming from TPCFastTransform (e.g. VDrift, T0..) are not set + template + static TPCFastTransformPOD& create(V& destVector, const TPCFastSpaceChargeCorrection& src); + + bool test(const TPCFastTransform& src, int npoints = 100000) const { return test(src.getCorrection(), npoints); } + bool test(const TPCFastSpaceChargeCorrection& origCorr, int npoints = 100000) const; +#endif + + static constexpr int NROWS = 152; + static constexpr int NSLICES = TPCFastTransformGeo::getNumberOfSlices(); + static constexpr int NSplineIDs = 3; ///< number of spline data sets for each slice/row + + private: +#if !defined(GPUCA_GPUCODE) + static constexpr size_t AlignmentBytes = 8; + static size_t alignOffset(size_t offs) + { + auto res = offs % AlignmentBytes; + return res ? offs + (AlignmentBytes - res) : offs; + } + static size_t estimateSize(const TPCFastTransform& src) { return estimateSize(src.getCorrection()); } + static size_t estimateSize(const TPCFastSpaceChargeCorrection& origCorr); + static TPCFastTransformPOD& create(char* buff, size_t buffSize, const TPCFastTransform& src); + static TPCFastTransformPOD& create(char* buff, size_t buffSize, const TPCFastSpaceChargeCorrection& src); + + ///< get address to which the offset in bytes must be added to arrive to particular dynamic part + GPUd() const char* getThis() const { return reinterpret_cast(this); } + GPUd() static TPCFastTransformPOD& getNonConst(char* head) { return *reinterpret_cast(head); } + +#endif + + /// Gives non-const TPC row info + GPUd() RowInfo& getRowInfo(int row) { return mRowInfo[row]; } + + ///< return offset of the spline object start (equivalent of mScenarioPtr in the TPCFastSpaceChargeCorrection) + GPUd() const size_t getScenarioOffset(int s) const { return (reinterpret_cast(getThis() + mOffsScenariosOffsets))[s]; } + + GPUd() void TransformInternal(int slice, int row, float& u, float& v, float& x) const; + + bool mApplyCorrections{}; ///< flag to apply corrections + int mNumberOfScenarios{}; ///< Number of approximation spline scenarios + size_t mTotalSize{}; ///< total size of the buffer + size_t mOffsScenariosOffsets{}; ///< start of the array of mNumberOfScenarios offsets for each type of spline + size_t mSplineDataOffsets[TPCFastTransformGeo::getNumberOfSlices()][NSplineIDs]; ///< start of data for each slice and iSpline data + long int mTimeStamp{}; ///< time stamp of the current calibration + float mInterpolationSafetyMargin{0.1f}; // 10% area around the TPC row. Outside of this area the interpolation returns the boundary values. + + /// _____ Parameters for drift length calculation ____ + /// + /// t = (float) time bin, y = global y + /// + /// L(t,y) = (t-mT0)*(mVdrift + mVdriftCorrY*y ) + mLdriftCorr ____ + /// + float mT0; ///< T0 in [time bin] + float mVdrift; ///< VDrift in [cm/time bin] + float mVdriftCorrY; ///< VDrift correction for global Y[cm] in [1/time bin] + float mLdriftCorr; ///< drift length correction in [cm] + + /// A coefficient for Time-Of-Flight correction: drift length -= EstimatedDistanceToVtx[cm]*mTOFcorr + /// + /// Since this correction requires a knowledge of the spatial position, it is appied after mCorrection, + /// not on the drift length but directly on V coordinate. + /// + /// mTOFcorr == mVdrift/(speed of light) + /// + float mTOFcorr; + + float mPrimVtxZ; ///< Z of the primary vertex, needed for the Time-Of-Flight correction + + TPCFastTransformGeo mGeo; ///< TPC geometry information + TPCFastSpaceChargeCorrection::SliceInfo mSliceInfo[TPCFastTransformGeo::getNumberOfSlices()]; ///< SliceInfo array + RowInfo mRowInfo[NROWS]; + SliceRowInfo mSliceRowInfo[NROWS * TPCFastTransformGeo::getNumberOfSlices()]; + + ClassDefNV(TPCFastTransformPOD, 0); +}; + +GPUdi() int TPCFastTransformPOD::getCorrection(int slice, int row, float u, float v, float& dx, float& du, float& dv) const +{ + const SplineType& spline = getSpline(slice, row); + const float* splineData = getSplineData(slice, row); + float gridU = 0, gridV = 0; + convUVtoGrid(slice, row, u, v, gridU, gridV); + float dxuv[3]; + spline.interpolateU(splineData, gridU, gridV, dxuv); + dx = dxuv[0]; + du = dxuv[1]; + dv = dxuv[2]; + return 0; +} + +GPUdi() int TPCFastTransformPOD::getCorrectionOld(int slice, int row, float u, float v, float& dx, float& du, float& dv) const +{ + const SplineType& spline = getSpline(slice, row); + const float* splineData = getSplineData(slice, row); + float gridU = 0, gridV = 0; + convUVtoGrid(slice, row, u, v, gridU, gridV); + float dxuv[3]; + spline.interpolateUold(splineData, gridU, gridV, dxuv); + dx = dxuv[0]; + du = dxuv[1]; + dv = dxuv[2]; + return 0; +} + +GPUdi() void TPCFastTransformPOD::getCorrectionInvCorrectedX(int slice, int row, float corrU, float corrV, float& x) const +{ + float gridU, gridV; + convCorrectedUVtoGrid(slice, row, corrU, corrV, gridU, gridV); + + const Spline2D& spline = reinterpret_cast&>(getSpline(slice, row)); + const float* splineData = getSplineData(slice, row, 1); + float dx = 0; + spline.interpolateU(splineData, gridU, gridV, &dx); + x = mGeo.getRowInfo(row).x + dx; +} + +GPUdi() void TPCFastTransformPOD::getCorrectionInvUV(int slice, int row, float corrU, float corrV, float& nomU, float& nomV) const +{ + float gridU, gridV; + convCorrectedUVtoGrid(slice, row, corrU, corrV, gridU, gridV); + + const Spline2D& spline = reinterpret_cast&>(getSpline(slice, row)); + const float* splineData = getSplineData(slice, row, 2); + + float duv[2]; + spline.interpolateU(splineData, gridU, gridV, duv); + nomU = corrU - duv[0]; + nomV = corrV - duv[1]; +} + +GPUdi() float TPCFastTransformPOD::getMaxDriftLength(int slice, int row, float pad) const +{ + const RowActiveArea& area = getSliceRowInfo(slice, row).activeArea; + const float* c = area.maxDriftLengthCheb; + float x = -1.f + 2.f * pad / mGeo.getRowInfo(row).maxPad; + float y = c[0] + c[1] * x; + float f0 = 1.f; + float f1 = x; + x *= 2.f; + for (int i = 2; i < 5; i++) { + double f = x * f1 - f0; + y += c[i] * f; + f0 = f1; + f1 = f; + } + return y; +} + +GPUdi() void TPCFastTransformPOD::schrinkUV(int slice, int row, float& u, float& v) const +{ + /// shrink u,v coordinats to the TPC row area +/- mInterpolationSafetyMargin + + float uWidth05 = mGeo.getRowInfo(row).getUwidth() * (0.5f + mInterpolationSafetyMargin); + float vWidth = mGeo.getTPCzLength(slice); + + if (u < -uWidth05) { + u = -uWidth05; + } + if (u > uWidth05) { + u = uWidth05; + } + if (v < -0.1f * vWidth) { + v = -0.1f * vWidth; + } + if (v > 1.1f * vWidth) { + v = 1.1f * vWidth; + } +} + +GPUdi() void TPCFastTransformPOD::schrinkCorrectedUV(int slice, int row, float& corrU, float& corrV) const +{ + /// shrink corrected u,v coordinats to the TPC row area +/- mInterpolationSafetyMargin + + const SliceRowInfo& sliceRowInfo = getSliceRowInfo(slice, row); + + float uMargin = mInterpolationSafetyMargin * mGeo.getRowInfo(row).getUwidth(); + float vMargin = mInterpolationSafetyMargin * mGeo.getTPCzLength(slice); + + if (corrU < sliceRowInfo.activeArea.cuMin - uMargin) { + corrU = sliceRowInfo.activeArea.cuMin - uMargin; + } + + if (corrU > sliceRowInfo.activeArea.cuMax + uMargin) { + corrU = sliceRowInfo.activeArea.cuMax + uMargin; + } + + if (corrV < 0.f - vMargin) { + corrV = 0.f - vMargin; + } + + if (corrV > sliceRowInfo.activeArea.cvMax + vMargin) { + corrV = sliceRowInfo.activeArea.cvMax + vMargin; + } +} + +GPUdi() void TPCFastTransformPOD::convUVtoGrid(int slice, int row, float u, float v, float& gu, float& gv) const +{ + // TODO optimise !!! + gu = 0.f; + gv = 0.f; + + schrinkUV(slice, row, u, v); + + const SliceRowInfo& info = getSliceRowInfo(slice, row); + const SplineType& spline = getSpline(slice, row); + + float su0 = 0.f, sv0 = 0.f; + mGeo.convUVtoScaledUV(slice, row, u, info.gridV0, su0, sv0); + mGeo.convUVtoScaledUV(slice, row, u, v, gu, gv); + + gv = (gv - sv0) / (1.f - sv0); + gu *= spline.getGridX1().getUmax(); + gv *= spline.getGridX2().getUmax(); +} + +GPUdi() void TPCFastTransformPOD::convGridToUV(int slice, int row, float gridU, float gridV, float& u, float& v) const +{ + // TODO optimise + /// convert u,v to internal grid coordinates + float su0 = 0.f, sv0 = 0.f; + const SliceRowInfo& info = getSliceRowInfo(slice, row); + const SplineType& spline = getSpline(slice, row); + mGeo.convUVtoScaledUV(slice, row, 0.f, info.gridV0, su0, sv0); + float su = gridU / spline.getGridX1().getUmax(); + float sv = sv0 + gridV / spline.getGridX2().getUmax() * (1.f - sv0); + mGeo.convScaledUVtoUV(slice, row, su, sv, u, v); +} + +GPUdi() void TPCFastTransformPOD::convCorrectedUVtoGrid(int slice, int row, float corrU, float corrV, float& gridU, float& gridV) const +{ + schrinkCorrectedUV(slice, row, corrU, corrV); + + const SliceRowInfo& sliceRowInfo = getSliceRowInfo(slice, row); + + gridU = (corrU - sliceRowInfo.gridCorrU0) * sliceRowInfo.scaleCorrUtoGrid; + gridV = (corrV - sliceRowInfo.gridCorrV0) * sliceRowInfo.scaleCorrVtoGrid; +} + +#if !defined(GPUCA_GPUCODE) +/// Create POD transform from old flat-buffer one. Provided vector will serve as a buffer +template +TPCFastTransformPOD& TPCFastTransformPOD::create(V& destVector, const TPCFastTransform& src) +{ + const auto& origCorr = src.getCorrection(); + size_t estSize = estimateSize(src); + destVector.resize(estSize); // allocate exact size + LOGP(debug, "OrigCorrSize:{} SelfSize: {} Estimated POS size: {}", src.getCorrection().getFlatBufferSize(), sizeof(TPCFastTransformPOD), estSize); + char* base = destVector.data(); + return create(destVector.data(), destVector.size(), src); +} + +template +TPCFastTransformPOD& TPCFastTransformPOD::create(V& destVector, const TPCFastSpaceChargeCorrection& origCorr) +{ + // create filling only part corresponding to TPCFastSpaceChargeCorrection. Data members coming from TPCFastTransform (e.g. VDrift, T0..) are not set + size_t estSize = estimateSize(origCorr); + destVector.resize(estSize); // allocate exact size + LOGP(debug, "OrigCorrSize:{} SelfSize: {} Estimated POS size: {}", origCorr.getFlatBufferSize(), sizeof(TPCFastTransformPOD), estSize); + char* base = destVector.data(); + return create(destVector.data(), destVector.size(), origCorr); +} +#endif + +GPUdi() void TPCFastTransformPOD::Transform(int slice, int row, float pad, float time, float& x, float& y, float& z, float vertexTime) const +{ + /// _______________ The main method: cluster transformation _______________________ + /// + /// Transforms raw TPC coordinates to local XYZ withing a slice + /// taking calibration + alignment into account. + const TPCFastTransformGeo::RowInfo& rowInfo = getGeometry().getRowInfo(row); + x = rowInfo.x; + float u = 0, v = 0; + convPadTimeToUV(slice, row, pad, time, u, v, vertexTime); + TransformInternal(slice, row, u, v, x); + getGeometry().convUVtoLocal(slice, u, v, y, z); + float dzTOF = 0; + getTOFcorrection(slice, row, x, y, z, dzTOF); + z += dzTOF; +} + +GPUdi() void TPCFastTransformPOD::TransformXYZ(int slice, int row, float& x, float& y, float& z) const +{ + float u, v; + getGeometry().convLocalToUV(slice, y, z, u, v); + TransformInternal(slice, row, u, v, x); + getGeometry().convUVtoLocal(slice, u, v, y, z); + float dzTOF = 0; + getTOFcorrection(slice, row, x, y, z, dzTOF); + z += dzTOF; +} + +GPUdi() void TPCFastTransformPOD::TransformInTimeFrame(int slice, float time, float& z, float maxTimeBin) const +{ + float v = 0; + convTimeToVinTimeFrame(slice, time, v, maxTimeBin); + getGeometry().convVtoLocal(slice, v, z); +} + +GPUdi() void TPCFastTransformPOD::TransformInTimeFrame(int slice, int row, float pad, float time, float& x, float& y, float& z, float maxTimeBin) const +{ + /// _______________ Special cluster transformation for a time frame _______________________ + /// + /// Same as Transform(), but clusters are shifted in z such, that Z(maxTimeBin)==0 + /// Corrections and Time-Of-Flight correction are not alpplied. + /// + const TPCFastTransformGeo::RowInfo& rowInfo = getGeometry().getRowInfo(row); + x = rowInfo.x; + float u = 0, v = 0; + convPadTimeToUVinTimeFrame(slice, row, pad, time, u, v, maxTimeBin); + getGeometry().convUVtoLocal(slice, u, v, y, z); +} + +GPUdi() void TPCFastTransformPOD::InverseTransformInTimeFrame(int slice, int row, float /*x*/, float y, float z, float& pad, float& time, float maxTimeBin) const +{ + /// Inverse transformation to TransformInTimeFrame + float u = 0, v = 0; + getGeometry().convLocalToUV(slice, y, z, u, v); + convUVtoPadTimeInTimeFrame(slice, row, u, v, pad, time, maxTimeBin); +} + +GPUdi() void TPCFastTransformPOD::InverseTransformYZtoX(int slice, int row, float y, float z, float& x) const +{ + /// Transformation y,z -> x + float u = 0, v = 0; + getGeometry().convLocalToUV(slice, y, z, u, v); + if (mApplyCorrections) { + getCorrectionInvCorrectedX(slice, row, u, v, x); + } else { + x = getGeometry().getRowInfo(row).x; // corrections are disabled + } + GPUCA_DEBUG_STREAMER_CHECK(if (o2::utils::DebugStreamer::checkStream(o2::utils::StreamFlags::streamFastTransform)) { + o2::utils::DebugStreamer::instance()->getStreamer("debug_fasttransform", "UPDATE") << o2::utils::DebugStreamer::instance()->getUniqueTreeName("tree_InverseTransformYZtoX").data() + << "slice=" << slice + << "row=" << row + << "y=" << y + << "z=" << z + << "x=" << x + << "v=" << v + << "u=" << u + << "\n"; + }) +} + +GPUdi() void TPCFastTransformPOD::InverseTransformYZtoNominalYZ(int slice, int row, float y, float z, float& ny, float& nz) const +{ + /// Transformation y,z -> x + float u = 0, v = 0, un = 0, vn = 0; + getGeometry().convLocalToUV(slice, y, z, u, v); + if (mApplyCorrections) { + getCorrectionInvUV(slice, row, u, v, un, vn); + } else { + un = u; + vn = v; + } + getGeometry().convUVtoLocal(slice, un, vn, ny, nz); + GPUCA_DEBUG_STREAMER_CHECK(if (o2::utils::DebugStreamer::checkStream(o2::utils::StreamFlags::streamFastTransform)) { + o2::utils::DebugStreamer::instance()->getStreamer("debug_fasttransform", "UPDATE") << o2::utils::DebugStreamer::instance()->getUniqueTreeName("tree_InverseTransformYZtoNominalYZ").data() + << "slice=" << slice + << "row=" << row + << "y=" << y + << "z=" << z + << "ny=" << ny + << "nz=" << nz + << "u=" << u + << "v=" << v + << "un=" << un + << "vn=" << vn + << "\n"; + }) +} + +GPUdi() void TPCFastTransformPOD::InverseTransformXYZtoNominalXYZ(int slice, int row, float x, float y, float z, float& nx, float& ny, float& nz) const +{ + /// Inverse transformation: Transformed X, Y and Z -> X, Y and Z, transformed w/o space charge correction + int row2 = row + 1; + if (row2 >= getGeometry().getNumberOfRows()) { + row2 = row - 1; + } + float nx1, ny1, nz1; // nominal coordinates for row + float nx2, ny2, nz2; // nominal coordinates for row2 + nx1 = getGeometry().getRowInfo(row).x; + nx2 = getGeometry().getRowInfo(row2).x; + InverseTransformYZtoNominalYZ(slice, row, y, z, ny1, nz1); + InverseTransformYZtoNominalYZ(slice, row2, y, z, ny2, nz2); + float c1 = (nx2 - nx) / (nx2 - nx1); + float c2 = (nx - nx1) / (nx2 - nx1); + nx = x; + ny = (ny1 * c1 + ny2 * c2); + nz = (nz1 * c1 + nz2 * c2); +} + +GPUdi() void TPCFastTransformPOD::TransformInternal(int slice, int row, float& u, float& v, float& x) const +{ + if (mApplyCorrections) { + float dx = 0.f, du = 0.f, dv = 0.f; + getCorrection(slice, row, u, v, dx, du, dv); + GPUCA_DEBUG_STREAMER_CHECK(if (o2::utils::DebugStreamer::checkStream(o2::utils::StreamFlags::streamFastTransform)) { + float ly, lz; + getGeometry().convUVtoLocal(slice, u, v, ly, lz); + float gx, gy, gz; + getGeometry().convLocalToGlobal(slice, x, ly, lz, gx, gy, gz); + float lyT, lzT; + float uCorr = u + du; + float vCorr = v + dv; + float lxT = x + dx; + getGeometry().convUVtoLocal(slice, uCorr, vCorr, lyT, lzT); + float invYZtoX; + InverseTransformYZtoX(slice, row, lyT, lzT, invYZtoX); + float YZtoNominalY; + float YZtoNominalZ; + InverseTransformYZtoNominalYZ(slice, row, lyT, lzT, YZtoNominalY, YZtoNominalZ); + o2::utils::DebugStreamer::instance()->getStreamer("debug_fasttransform", "UPDATE") << o2::utils::DebugStreamer::instance()->getUniqueTreeName("tree_Transform").data() + // corrections in x, u, v + << "dx=" << dx + << "du=" << du + << "dv=" << dv + << "v=" << v + << "u=" << u + << "row=" << row + << "slice=" << slice + // original local coordinates + << "ly=" << ly + << "lz=" << lz + << "lx=" << x + // corrected local coordinated + << "lxT=" << lxT + << "lyT=" << lyT + << "lzT=" << lzT + // global uncorrected coordinates + << "gx=" << gx + << "gy=" << gy + << "gz=" << gz + // some transformations which are applied + << "invYZtoX=" << invYZtoX + << "YZtoNominalY=" << YZtoNominalY + << "YZtoNominalZ=" << YZtoNominalZ + << "\n"; + }) + x += dx; + u += du; + v += dv; + } +} + +GPUdi() void TPCFastTransformPOD::TransformIdealZ(int slice, float time, float& z, float vertexTime) const +{ + /// _______________ The main method: cluster transformation _______________________ + /// + /// Transforms time TPC coordinates to local Z withing a slice + /// Ideal transformation: only Vdrift from DCS. + /// No space charge corrections, no time of flight correction + /// + + float v = (time - mT0 - vertexTime) * mVdrift; // drift length cm + getGeometry().convVtoLocal(slice, v, z); +} + +GPUdi() void TPCFastTransformPOD::TransformIdeal(int slice, int row, float pad, float time, float& x, float& y, float& z, float vertexTime) const +{ + /// _______________ The main method: cluster transformation _______________________ + /// + /// Transforms raw TPC coordinates to local XYZ withing a slice + /// Ideal transformation: only Vdrift from DCS. + /// No space charge corrections, no time of flight correction + /// + const TPCFastTransformGeo::RowInfo& rowInfo = getGeometry().getRowInfo(row); + x = rowInfo.x; + float u = (pad - 0.5f * rowInfo.maxPad) * rowInfo.padWidth; + float v = (time - mT0 - vertexTime) * mVdrift; // drift length cm + getGeometry().convUVtoLocal(slice, u, v, y, z); +} + +GPUdi() void TPCFastTransformPOD::convPadTimeToUV(int slice, int row, float pad, float time, float& u, float& v, float vertexTime) const +{ + bool sideC = (slice >= getGeometry().getNumberOfSlicesA()); + const TPCFastTransformGeo::RowInfo& rowInfo = getGeometry().getRowInfo(row); + const TPCFastTransformGeo::SliceInfo& sliceInfo = getGeometry().getSliceInfo(slice); + float x = rowInfo.x; + u = (pad - 0.5f * rowInfo.maxPad) * rowInfo.padWidth; + float y = sideC ? -u : u; // pads are mirrorred on C-side + float yLab = y * sliceInfo.cosAlpha + x * sliceInfo.sinAlpha; + v = (time - mT0 - vertexTime) * (mVdrift + mVdriftCorrY * yLab) + mLdriftCorr; // drift length cm +} + +GPUdi() void TPCFastTransformPOD::convTimeToVinTimeFrame(int slice, float time, float& v, float maxTimeBin) const +{ + v = (time - mT0 - maxTimeBin) * mVdrift + mLdriftCorr; // drift length cm + if (slice < getGeometry().getNumberOfSlicesA()) { + v += getGeometry().getTPCzLengthA(); + } else { + v += getGeometry().getTPCzLengthC(); + } +} + +GPUdi() void TPCFastTransformPOD::convPadTimeToUVinTimeFrame(int slice, int row, float pad, float time, float& u, float& v, float maxTimeBin) const +{ + const TPCFastTransformGeo::RowInfo& rowInfo = getGeometry().getRowInfo(row); + u = (pad - 0.5f * rowInfo.maxPad) * rowInfo.padWidth; + convTimeToVinTimeFrame(slice, time, v, maxTimeBin); +} + +GPUdi() float TPCFastTransformPOD::convZOffsetToVertexTime(int slice, float zOffset, float maxTimeBin) const +{ + if (slice < getGeometry().getNumberOfSlicesA()) { + return maxTimeBin - (getGeometry().getTPCzLengthA() + zOffset) / mVdrift; + } else { + return maxTimeBin - (getGeometry().getTPCzLengthC() - zOffset) / mVdrift; + } +} + +GPUdi() float TPCFastTransformPOD::convVertexTimeToZOffset(int slice, float vertexTime, float maxTimeBin) const +{ + if (slice < getGeometry().getNumberOfSlicesA()) { + return (maxTimeBin - vertexTime) * mVdrift - getGeometry().getTPCzLengthA(); + } else { + return -((maxTimeBin - vertexTime) * mVdrift - getGeometry().getTPCzLengthC()); + } +} + +GPUdi() void TPCFastTransformPOD::convUVtoPadTime(int slice, int row, float u, float v, float& pad, float& time, float vertexTime) const +{ + bool sideC = (slice >= getGeometry().getNumberOfSlicesA()); + + const TPCFastTransformGeo::RowInfo& rowInfo = getGeometry().getRowInfo(row); + const TPCFastTransformGeo::SliceInfo& sliceInfo = getGeometry().getSliceInfo(slice); + + pad = u / rowInfo.padWidth + 0.5f * rowInfo.maxPad; + + float x = rowInfo.x; + float y = sideC ? -u : u; // pads are mirrorred on C-side + float yLab = y * sliceInfo.cosAlpha + x * sliceInfo.sinAlpha; + time = mT0 + vertexTime + (v - mLdriftCorr) / (mVdrift + mVdriftCorrY * yLab); +} + +GPUdi() void TPCFastTransformPOD::convVtoTime(float v, float& time, float vertexTime) const +{ + float yLab = 0.f; + time = mT0 + vertexTime + (v - mLdriftCorr) / (mVdrift + mVdriftCorrY * yLab); +} + +GPUdi() void TPCFastTransformPOD::convUVtoPadTimeInTimeFrame(int slice, int row, float u, float v, float& pad, float& time, float maxTimeBin) const +{ + if (slice < getGeometry().getNumberOfSlicesA()) { + v -= getGeometry().getTPCzLengthA(); + } else { + v -= getGeometry().getTPCzLengthC(); + } + const TPCFastTransformGeo::RowInfo& rowInfo = getGeometry().getRowInfo(row); + pad = u / rowInfo.padWidth + 0.5f * rowInfo.maxPad; + time = mT0 + maxTimeBin + (v - mLdriftCorr) / mVdrift; +} + +GPUdi() void TPCFastTransformPOD::getTOFcorrection(int slice, int /*row*/, float x, float y, float z, float& dz) const +{ + // calculate time of flight correction for z coordinate + + bool sideC = (slice >= getGeometry().getNumberOfSlicesA()); + float distZ = z - mPrimVtxZ; + float dv = -GPUCommonMath::Sqrt(x * x + y * y + distZ * distZ) * mTOFcorr; + dz = sideC ? dv : -dv; +} + +GPUdi() float TPCFastTransformPOD::getMaxDriftTime(int slice, int row, float pad) const +{ + /// maximal possible drift time of the active area + float maxL = getMaxDriftLength(slice, row, pad); + bool sideC = (slice >= getGeometry().getNumberOfSlicesA()); + const TPCFastTransformGeo::RowInfo& rowInfo = getGeometry().getRowInfo(row); + const TPCFastTransformGeo::SliceInfo& sliceInfo = getGeometry().getSliceInfo(slice); + float x = rowInfo.x; + float u = (pad - 0.5f * rowInfo.maxPad) * rowInfo.padWidth; + + float y = sideC ? -u : u; // pads are mirrorred on C-side + float yLab = y * sliceInfo.cosAlpha + x * sliceInfo.sinAlpha; + return mT0 + (maxL - mLdriftCorr) / (mVdrift + mVdriftCorrY * yLab); +} + +GPUdi() float TPCFastTransformPOD::getMaxDriftTime(int slice, int row) const +{ + /// maximal possible drift time of the active area + float maxL = getMaxDriftLength(slice, row); + float maxTime = 0.f; + convVtoTime(maxL, maxTime, 0.f); + return maxTime; +} + +GPUdi() float TPCFastTransformPOD::getMaxDriftTime(int slice) const +{ + /// maximal possible drift time of the active area + float maxL = getMaxDriftLength(slice); + float maxTime = 0.f; + convVtoTime(maxL, maxTime, 0.f); + return maxTime; +} + +} // namespace gpu +} // namespace GPUCA_NAMESPACE + +#endif // ALICEO2_GPU_TPCFastTransformPOD_H diff --git a/GPU/TPCFastTransformation/TPCFastTransformationLinkDef_O2.h b/GPU/TPCFastTransformation/TPCFastTransformationLinkDef_O2.h index 4421d44aab0c8..7482cbf52768a 100644 --- a/GPU/TPCFastTransformation/TPCFastTransformationLinkDef_O2.h +++ b/GPU/TPCFastTransformation/TPCFastTransformationLinkDef_O2.h @@ -63,6 +63,7 @@ #pragma link C++ class o2::gpu::TPCFastTransformGeo::RowInfo + ; #pragma link C++ class o2::gpu::TPCFastTransform + ; +#pragma link C++ class o2::gpu::TPCFastTransformPOD; #pragma link C++ class o2::gpu::TPCFastSpaceChargeCorrectionMap + ; #pragma link C++ class o2::gpu::TPCFastSpaceChargeCorrection::RowInfo + ; From 73187e41f94b2c7c5d56fc6243234caceff5c012 Mon Sep 17 00:00:00 2001 From: shahoian Date: Thu, 18 Jul 2024 15:56:54 +0200 Subject: [PATCH 2/2] Add LumiIDC data member to TPCFastTransform --- GPU/TPCFastTransformation/TPCFastTransform.cxx | 5 ++++- GPU/TPCFastTransformation/TPCFastTransform.h | 7 ++++++- 2 files changed, 10 insertions(+), 2 deletions(-) diff --git a/GPU/TPCFastTransformation/TPCFastTransform.cxx b/GPU/TPCFastTransformation/TPCFastTransform.cxx index 6424292a6035a..855c63cc2fe2e 100644 --- a/GPU/TPCFastTransformation/TPCFastTransform.cxx +++ b/GPU/TPCFastTransformation/TPCFastTransform.cxx @@ -37,7 +37,7 @@ using namespace GPUCA_NAMESPACE::gpu; TPCFastTransform::TPCFastTransform() - : FlatObject(), mTimeStamp(0), mCorrection(), mApplyCorrection(1), mT0(0.f), mVdrift(0.f), mVdriftCorrY(0.f), mLdriftCorr(0.f), mTOFcorr(0.f), mPrimVtxZ(0.f), mLumi(0.f), mLumiError(0.f), mLumiScaleFactor(1.0f) + : FlatObject(), mTimeStamp(0), mCorrection(), mApplyCorrection(1), mT0(0.f), mVdrift(0.f), mVdriftCorrY(0.f), mLdriftCorr(0.f), mTOFcorr(0.f), mPrimVtxZ(0.f), mLumiIDC(0.f), mLumi(0.f), mLumiError(0.f), mLumiScaleFactor(1.0f) { // Default Constructor: creates an empty uninitialized object } @@ -58,6 +58,7 @@ void TPCFastTransform::cloneFromObject(const TPCFastTransform& obj, char* newFla mLdriftCorr = obj.mLdriftCorr; mTOFcorr = obj.mTOFcorr; mPrimVtxZ = obj.mPrimVtxZ; + mLumiIDC = obj.mLumiIDC; mLumi = obj.mLumi; mLumiError = obj.mLumiError; mLumiScaleFactor = obj.mLumiScaleFactor; @@ -108,6 +109,7 @@ void TPCFastTransform::startConstruction(const TPCFastSpaceChargeCorrection& cor mLdriftCorr = 0.f; mTOFcorr = 0.f; mPrimVtxZ = 0.f; + mLumiIDC = 0.f; mLumi = 0.f; mLumiError = 0.f; mLumiScaleFactor = 1.f; @@ -158,6 +160,7 @@ void TPCFastTransform::print() const LOG(info) << "mLdriftCorr = " << mLdriftCorr; LOG(info) << "mTOFcorr = " << mTOFcorr; LOG(info) << "mPrimVtxZ = " << mPrimVtxZ; + LOG(info) << "mLumiIDC = " << mLumiIDC; LOG(info) << "mLumi = " << mLumi; LOG(info) << "mLumiError = " << mLumiError; LOG(info) << "mLumiScaleFactor = " << mLumiScaleFactor; diff --git a/GPU/TPCFastTransformation/TPCFastTransform.h b/GPU/TPCFastTransformation/TPCFastTransform.h index 492702a185c08..c431f8459367c 100644 --- a/GPU/TPCFastTransformation/TPCFastTransform.h +++ b/GPU/TPCFastTransformation/TPCFastTransform.h @@ -162,6 +162,7 @@ class TPCFastTransform : public FlatObject /// Set Lumi info void setLumi(float l) { mLumi = l; } + void setLumiIDC(float l) { mLumiIDC = l; } void setLumiError(float e) { mLumiError = e; } void setLumiScaleFactor(float s) { mLumiScaleFactor = s; } @@ -256,6 +257,9 @@ class TPCFastTransform : public FlatObject /// Return map lumi GPUd() float getLumi() const { return mLumi; } + /// Return map lumiIDC + GPUd() float getLumiIDC() const { return mLumiIDC; } + /// Return map lumi error GPUd() float getLumiError() const { return mLumiError; } @@ -335,6 +339,7 @@ class TPCFastTransform : public FlatObject float mPrimVtxZ; ///< Z of the primary vertex, needed for the Time-Of-Flight correction + float mLumiIDC; ///< luminosity estimator in IDC units float mLumi; ///< luminosity estimator float mLumiError; ///< error on luminosity float mLumiScaleFactor; ///< user correction factor for lumi (e.g. normalization, efficiency correction etc.) @@ -345,7 +350,7 @@ class TPCFastTransform : public FlatObject GPUd() void TransformInternal(int slice, int row, float& u, float& v, float& x, const TPCFastTransform* ref, const TPCFastTransform* ref2, float scale, float scale2, int scaleMode) const; #ifndef GPUCA_ALIROOT_LIB - ClassDefNV(TPCFastTransform, 3); + ClassDefNV(TPCFastTransform, 4); #endif };