Skip to content

Commit

Permalink
TPC Splines: add limits for SP correction values per TPC row
Browse files Browse the repository at this point in the history
  • Loading branch information
@sgorbuno
sgorbuno committed Jun 20, 2024
1 parent fa1d8d7 commit dec1e31
Show file tree
Hide file tree
Showing 3 changed files with 97 additions and 21 deletions.
18 changes: 17 additions & 1 deletion Detectors/TPC/calibration/src/TPCFastSpaceChargeCorrectionHelper.cxx
View file
Original file line number Diff line number Diff line change
Expand Up @@ -144,6 +144,9 @@ void TPCFastSpaceChargeCorrectionHelper::fillSpaceChargeCorrectionFromMap(TPCFas
float* splineParameters = correction.getSplineData(slice, row);
const std::vector<o2::gpu::TPCFastSpaceChargeCorrectionMap::CorrectionPoint>& data = mCorrectionMap.getPoints(slice, row);
int nDataPoints = data.size();
auto& info = correction.getSliceRowInfo(slice, row);
info.resetMaxValues();
info.resetMaxValuesInv();
if (nDataPoints >= 4) {
std::vector<double> pointSU(nDataPoints);
std::vector<double> pointSV(nDataPoints);
Expand All @@ -156,6 +159,8 @@ void TPCFastSpaceChargeCorrectionHelper::fillSpaceChargeCorrectionFromMap(TPCFas
pointCorr[3 * i + 0] = dx;
pointCorr[3 * i + 1] = du;
pointCorr[3 * i + 2] = dv;
info.updateMaxValues(2. * dx, 2. * du, 2. * dv);
info.updateMaxValuesInv(-2. * dx, -2. * du, -2. * dv);
}
helper.approximateDataPoints(spline, splineParameters, 0., spline.getGridX1().getUmax(), 0., spline.getGridX2().getUmax(), &pointSU[0],
&pointSV[0], &pointCorr[0], nDataPoints);
Expand Down Expand Up @@ -767,9 +772,20 @@ std::unique_ptr<o2::gpu::TPCFastSpaceChargeCorrection> TPCFastSpaceChargeCorrect

double yStep = (yLast - yFirst) / 2;

double zFirst = z - dz / 2.;
double zLast = z + dz / 2.;
double zStep = (zLast - zFirst) / 2.;

if (0) { // no smoothing
yFirst = y;
yLast = y;
zFirst = z;
zLast = z;
}

for (double py = yFirst; py <= yLast + yStep / 2.; py += yStep) {

for (double pz = z - dz / 2.; pz <= z + dz / 2. + 1.e-4; pz += dz / 2.) {
for (double pz = zFirst; pz <= zLast + zStep / 2.; pz += zStep) {
map.addCorrectionPoint(iRoc, iRow, py, pz, correctionX, correctionY,
correctionZ);
}
Expand Down
85 changes: 71 additions & 14 deletions GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,6 +21,7 @@
#include "TPCFastTransformGeo.h"
#include "FlatObject.h"
#include "GPUCommonDef.h"
#include "GPUCommonMath.h"

namespace GPUCA_NAMESPACE
{
Expand Down Expand Up @@ -61,15 +62,64 @@ class TPCFastSpaceChargeCorrection : public FlatObject
};

struct SliceRowInfo {
float gridU0{0.f}; //< U coordinate of the U-grid start
float scaleUtoGrid{0.f}; //< scale U to U-grid coordinate
float gridV0{0.f}; ///< V coordinate of the V-grid start
float scaleVtoGrid{0.f}; //< scale V to V-grid coordinate
float gridCorrU0{0.f}; ///< U coordinate of the U-grid start for corrected U
float scaleCorrUtoGrid{0.f}; ///< scale corrected U to U-grid coordinate
float gridCorrV0{0.f}; ///< V coordinate of the V-grid start for corrected V
float scaleCorrVtoGrid{0.f}; ///< scale corrected V to V-grid coordinate
float gridU0{0.f}; //< U coordinate of the U-grid start
float scaleUtoGrid{0.f}; //< scale U to U-grid coordinate
float gridV0{0.f}; ///< V coordinate of the V-grid start
float scaleVtoGrid{0.f}; //< scale V to V-grid coordinate
float gridCorrU0{0.f}; ///< U coordinate of the U-grid start for corrected U
float scaleCorrUtoGrid{0.f}; ///< scale corrected U to U-grid coordinate
float gridCorrV0{0.f}; ///< V coordinate of the V-grid start for corrected V
float scaleCorrVtoGrid{0.f}; ///< scale corrected V to V-grid coordinate
float maxCorr[3]{10.f, 10.f, 10.f}; ///< max correction for dX, dU, dV
float minCorr[3]{-10.f, -10.f, -10.f}; ///< min correction for dX, dU, dV
float maxInvCorr[3]{10.f, 10.f, 10.f}; ///< max inverse correction for dX, dU, dV
float minInvCorr[3]{-10.f, -10.f, -10.f}; ///< min inverse correction for dX, dU, dV
RowActiveArea activeArea;

void resetMaxValues()
{
maxCorr[0] = 1.f;
minCorr[0] = -1.f;
maxCorr[1] = 1.f;
minCorr[1] = -1.f;
maxCorr[2] = 1.f;
minCorr[2] = -1.f;
}

void updateMaxValues(float dx, float du, float dv)
{
maxCorr[0] = GPUCommonMath::Max(maxCorr[0], dx);
minCorr[0] = GPUCommonMath::Min(minCorr[0], dx);

maxCorr[1] = GPUCommonMath::Max(maxCorr[1], du);
minCorr[1] = GPUCommonMath::Min(minCorr[1], du);

maxCorr[2] = GPUCommonMath::Max(maxCorr[2], dv);
minCorr[2] = GPUCommonMath::Min(minCorr[2], dv);
}

void resetMaxValuesInv()
{
maxInvCorr[0] = 1.f;
minInvCorr[0] = -1.f;
maxInvCorr[1] = 1.f;
minInvCorr[1] = -1.f;
maxInvCorr[2] = 1.f;
minInvCorr[2] = -1.f;
}

void updateMaxValuesInv(float dx, float du, float dv)
{
maxInvCorr[0] = GPUCommonMath::Max(maxInvCorr[0], dx);
minInvCorr[0] = GPUCommonMath::Min(minInvCorr[0], dx);

maxInvCorr[1] = GPUCommonMath::Max(maxInvCorr[1], du);
minInvCorr[1] = GPUCommonMath::Min(minInvCorr[1], du);

maxInvCorr[2] = GPUCommonMath::Max(maxInvCorr[2], dv);
minInvCorr[2] = GPUCommonMath::Min(minInvCorr[2], dv);
}

#ifndef GPUCA_ALIROOT_LIB
ClassDefNV(SliceRowInfo, 2);
#endif
Expand Down Expand Up @@ -406,9 +456,10 @@ GPUdi() int TPCFastSpaceChargeCorrection::getCorrection(int slice, int row, floa
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];
const auto& info = getSliceRowInfo(slice, row);
dx = GPUCommonMath::Max(info.minCorr[0], GPUCommonMath::Min(info.maxCorr[0], dxuv[0]));
du = GPUCommonMath::Max(info.minCorr[1], GPUCommonMath::Min(info.maxCorr[1], dxuv[1]));
dv = GPUCommonMath::Max(info.minCorr[2], GPUCommonMath::Min(info.maxCorr[2], dxuv[2]));
return 0;
}

Expand All @@ -420,9 +471,10 @@ GPUdi() int TPCFastSpaceChargeCorrection::getCorrectionOld(int slice, int row, f
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];
const auto& info = getSliceRowInfo(slice, row);
dx = GPUCommonMath::Max(info.minCorr[0], GPUCommonMath::Min(info.maxCorr[0], dxuv[0]));
du = GPUCommonMath::Max(info.minCorr[1], GPUCommonMath::Min(info.maxCorr[1], dxuv[1]));
dv = GPUCommonMath::Max(info.minCorr[2], GPUCommonMath::Min(info.maxCorr[2], dxuv[2]));
return 0;
}

Expand All @@ -436,6 +488,8 @@ GPUdi() void TPCFastSpaceChargeCorrection::getCorrectionInvCorrectedX(
const float* splineData = getSplineData(slice, row, 1);
float dx = 0;
spline.interpolateU(splineData, gridU, gridV, &dx);
const auto& info = getSliceRowInfo(slice, row);
dx = GPUCommonMath::Max(info.minInvCorr[0], GPUCommonMath::Min(info.maxInvCorr[0], dx));
x = mGeo.getRowInfo(row).x + dx;
}

Expand All @@ -450,6 +504,9 @@ GPUdi() void TPCFastSpaceChargeCorrection::getCorrectionInvUV(

float duv[2];
spline.interpolateU(splineData, gridU, gridV, duv);
const auto& info = getSliceRowInfo(slice, row);
duv[0] = GPUCommonMath::Max(info.minInvCorr[1], GPUCommonMath::Min(info.maxInvCorr[1], duv[0]));
duv[1] = GPUCommonMath::Max(info.minInvCorr[2], GPUCommonMath::Min(info.maxInvCorr[2], duv[1]));
nomU = corrU - duv[0];
nomV = corrV - duv[1];
}
Expand Down
15 changes: 9 additions & 6 deletions GPU/TPCFastTransformation/macro/TPCFastTransformInit.C
View file
Original file line number Diff line number Diff line change
Expand Up @@ -21,8 +21,6 @@
/// root -l TPCFastTransformInit.C'("debugVoxRes.root")'
///

#include "Algorithm/RangeTokenizer.h"

#if !defined(__CLING__) || defined(__ROOTCLING__)

#include <filesystem>
Expand All @@ -41,6 +39,8 @@
#include "TPCCalibration/TPCFastSpaceChargeCorrectionHelper.h"
#endif

#include "Algorithm/RangeTokenizer.h"

using namespace o2::tpc;
using namespace o2::gpu;

Expand Down Expand Up @@ -99,8 +99,9 @@ void TPCFastTransformInit(const char* fileName = "debugVoxRes.root",
trackResiduals.setZ2XBinning(z2xBins);
trackResiduals.init();

{
std::cout << "input track residuals: " << std::endl;
{ // debug output

std::cout << " ===== input track residuals ==== " << std::endl;
std::cout << "voxel tree y2xBins: " << y2xBins.size() << std::endl;

for (auto y2x : y2xBins) {
Expand All @@ -127,6 +128,7 @@ void TPCFastTransformInit(const char* fileName = "debugVoxRes.root",
for (int i = 0; i < nZ2Xbins; i++) {
std::cout << "getZ2X(bin) : " << trackResiduals.getZ2X(i) << std::endl;
}
std::cout << " ==================================== " << std::endl;
}

std::cout << "create fast transformation ... " << std::endl;
Expand Down Expand Up @@ -310,15 +312,16 @@ void TPCFastTransformInit(const char* fileName = "debugVoxRes.root",
geo.convUVtoLocal(iRoc, u + cu, v + cv, cy, cz);
cy -= y;
cz -= z;

double d[3] = {cx - correctionX, cy - correctionY, cz - correctionZ};
if (voxEntries >= 1.) {
for (int i = 0; i < 3; i++) {
if (fabs(maxDiff[i]) < fabs(d[i])) {
maxDiff[i] = d[i];
maxDiffRoc[i] = iRoc;
maxDiffRow[i] = iRow;
std::cout << " roc " << iRoc << " row " << iRow << " xyz " << i
<< " diff " << d[i] << " entries " << voxEntries << " y " << y2xBin << " z " << z2xBin << std::endl;
// std::cout << " roc " << iRoc << " row " << iRow << " xyz " << i
// << " diff " << d[i] << " entries " << voxEntries << " y " << y2xBin << " z " << z2xBin << std::endl;
}
sumDiff[i] += d[i] * d[i];
}
Expand Down

0 comments on commit dec1e31

Please sign in to comment.