diff --git a/Detectors/TPC/calibration/src/TPCFastSpaceChargeCorrectionHelper.cxx b/Detectors/TPC/calibration/src/TPCFastSpaceChargeCorrectionHelper.cxx index 783c1837590b9..5a26dabaa2db5 100644 --- a/Detectors/TPC/calibration/src/TPCFastSpaceChargeCorrectionHelper.cxx +++ b/Detectors/TPC/calibration/src/TPCFastSpaceChargeCorrectionHelper.cxx @@ -137,19 +137,13 @@ void TPCFastSpaceChargeCorrectionHelper::fillSpaceChargeCorrectionFromMap(TPCFas auto myThread = [&](int iThread) { for (int row = iThread; row < correction.getGeometry().getNumberOfRows(); row += mNthreads) { - TPCFastSpaceChargeCorrection::SplineType& spline = correction.getSpline(sector, row); + TPCFastSpaceChargeCorrection::SplineType& spline = correction.getSplineForRow(row); Spline2DHelper helper; std::vector splineParameters; splineParameters.resize(spline.getNumberOfParameters()); const std::vector& data = mCorrectionMap.getPoints(sector, row); int nDataPoints = data.size(); - auto& info = correction.getSectorRowInfo(sector, row); - if (!processingInverseCorrection) { - info.resetMaxValues(); - } - info.updateMaxValues(1., 1., 1.); - info.updateMaxValues(-1., -1., -1.); if (nDataPoints >= 4) { std::vector pointGU(nDataPoints); @@ -170,7 +164,6 @@ void TPCFastSpaceChargeCorrectionHelper::fillSpaceChargeCorrectionFromMap(TPCFas pointCorr[3 * i + 0] = p.mDx; pointCorr[3 * i + 1] = p.mDy; pointCorr[3 * i + 2] = p.mDz; - info.updateMaxValues(5. * p.mDx, 5. * p.mDy, 5. * p.mDz); } helper.approximateDataPoints(spline, splineParameters.data(), 0., spline.getGridX1().getUmax(), 0., spline.getGridX2().getUmax(), pointGU.data(), pointGV.data(), pointCorr.data(), pointWeight.data(), nDataPoints); @@ -263,14 +256,12 @@ std::unique_ptr TPCFastSpaceChargeCorrectionHelper correction.startConstruction(mGeo, nCorrectionScenarios); // assign spline type for TPC rows - for (int sector = 0; sector < mGeo.getNumberOfSectors(); sector++) { - for (int row = 0; row < mGeo.getNumberOfRows(); row++) { - int scenario = row / 10; - if (scenario >= nCorrectionScenarios) { - scenario = nCorrectionScenarios - 1; - } - correction.setRowScenarioID(sector, row, scenario); + for (int row = 0; row < mGeo.getNumberOfRows(); row++) { + int scenario = row / 10; + if (scenario >= nCorrectionScenarios) { + scenario = nCorrectionScenarios - 1; } + correction.setRowScenarioID(row, scenario); } for (int scenario = 0; scenario < nCorrectionScenarios; scenario++) { @@ -408,11 +399,10 @@ std::unique_ptr TPCFastSpaceChargeCorrect int nY2Xbins = trackResiduals.getNY2XBins(); int nZ2Xbins = trackResiduals.getNZ2XBins(); - std::vector knotsDouble[3]; + std::vector knotsDouble[2]; knotsDouble[0].reserve(nY2Xbins); knotsDouble[1].reserve(nZ2Xbins); - knotsDouble[2].reserve(nZ2Xbins); // to get enouth measurements, make a spline knot at every second bin. Boundary bins are always included. @@ -425,16 +415,14 @@ std::unique_ptr TPCFastSpaceChargeCorrect for (int i = 0, j = nZ2Xbins - 1; i <= j; i += 2, j -= 2) { knotsDouble[1].push_back(trackResiduals.getZ2X(i)); - knotsDouble[2].push_back(-trackResiduals.getZ2X(i)); if (j >= i + 1) { knotsDouble[1].push_back(trackResiduals.getZ2X(j)); - knotsDouble[2].push_back(-trackResiduals.getZ2X(j)); } } - std::vector knotsInt[3]; + std::vector knotsInt[2]; - for (int dim = 0; dim < 3; dim++) { + for (int dim = 0; dim < 2; dim++) { auto& knotsD = knotsDouble[dim]; std::sort(knotsD.begin(), knotsD.end()); @@ -470,12 +458,10 @@ std::unique_ptr TPCFastSpaceChargeCorrect } auto& yKnotsInt = knotsInt[0]; - auto& zKnotsIntA = knotsInt[1]; - auto& zKnotsIntC = knotsInt[2]; + auto& zKnotsInt = knotsInt[1]; int nKnotsY = yKnotsInt.size(); - int nKnotsZA = zKnotsIntA.size(); - int nKnotsZC = zKnotsIntC.size(); + int nKnotsZ = zKnotsInt.size(); // std::cout << "n knots Y: " << nKnotsY << std::endl; // std::cout << "n knots Z: " << nKnotsZA << ", " << nKnotsZC << std::endl; @@ -485,53 +471,41 @@ std::unique_ptr TPCFastSpaceChargeCorrect { // create the correction object - const int nCorrectionScenarios = 2; // different grids for TPC A and TPC C sides + const int nCorrectionScenarios = 1; correction.startConstruction(geo, nCorrectionScenarios); // init rows - for (int iSector = 0; iSector < nSectors; iSector++) { - int id = iSector < geo.getNumberOfSectorsA() ? 0 : 1; - for (int row = 0; row < geo.getNumberOfRows(); row++) { - correction.setRowScenarioID(iSector, row, id); - } + for (int row = 0; row < geo.getNumberOfRows(); row++) { + correction.setRowScenarioID(row, 0); } + { // init spline scenario TPCFastSpaceChargeCorrection::SplineType spline; - spline.recreate(nKnotsY, &yKnotsInt[0], nKnotsZA, &zKnotsIntA[0]); + spline.recreate(nKnotsY, &yKnotsInt[0], nKnotsZ, &zKnotsInt[0]); correction.setSplineScenario(0, spline); - spline.recreate(nKnotsY, &yKnotsInt[0], nKnotsZC, &zKnotsIntC[0]); - correction.setSplineScenario(1, spline); } correction.finishConstruction(); } // .. create the correction object // set the grid borders - for (int iSector = 0; iSector < geo.getNumberOfSectors(); iSector++) { - for (int iRow = 0; iRow < geo.getNumberOfRows(); iRow++) { - auto& info = correction.getSectorRowInfo(iSector, iRow); - const auto& spline = correction.getSpline(iSector, iRow); - double rowX = geo.getRowInfo(iRow).x; - double yMin = rowX * trackResiduals.getY2X(iRow, 0); - double yMax = rowX * trackResiduals.getY2X(iRow, trackResiduals.getNY2XBins() - 1); - double zMin = rowX * trackResiduals.getZ2X(0); - double zMax = rowX * trackResiduals.getZ2X(trackResiduals.getNZ2XBins() - 1); - double zOut = zMax; - if (iSector >= geo.getNumberOfSectorsA()) { - // TPC C side - zOut = -zOut; - zMax = -zMin; - zMin = zOut; - } - info.gridMeasured.set(yMin, spline.getGridX1().getUmax() / (yMax - yMin), // y - zMin, spline.getGridX2().getUmax() / (zMax - zMin), // z - zOut, geo.getZreadout(iSector)); // correction scaling region - - info.gridReal = info.gridMeasured; - - // std::cout << " iSector " << iSector << " iRow " << iRow << " uMin: " << uMin << " uMax: " << uMax << " vMin: " << vMin << " vMax: " << vMax - //<< " grid scale u "<< info.scaleUtoGrid << " grid scale v "<< info.scaleVtoGrid<< std::endl; - } + for (int iRow = 0; iRow < geo.getNumberOfRows(); iRow++) { + auto& info = correction.getRowInfo(iRow); + const auto& spline = correction.getSplineForRow(iRow); + double rowX = geo.getRowInfo(iRow).x; + double yMin = rowX * trackResiduals.getY2X(iRow, 0); + double yMax = rowX * trackResiduals.getY2X(iRow, trackResiduals.getNY2XBins() - 1); + double zMin = rowX * trackResiduals.getZ2X(0); + double zMax = rowX * trackResiduals.getZ2X(trackResiduals.getNZ2XBins() - 1); + double zOut = zMax; + info.gridMeasured.set(yMin, spline.getGridX1().getUmax() / (yMax - yMin), // y + zMin, spline.getGridX2().getUmax() / (zMax - zMin), // z + zOut, geo.getTPCzLength()); // correction scaling region + + info.gridReal = info.gridMeasured; + + // std::cout << " iSector " << iSector << " iRow " << iRow << " uMin: " << uMin << " uMax: " << uMax << " vMin: " << vMin << " vMax: " << vMax + //<< " grid scale u "<< info.scaleUtoGrid << " grid scale v "<< info.scaleVtoGrid<< std::endl; } LOG(info) << "fast space charge correction helper: preparation took " << watch1.RealTime() << "s"; @@ -783,8 +757,8 @@ std::unique_ptr TPCFastSpaceChargeCorrect // feed the row data to the helper - auto& info = correction.getSectorRowInfo(iSector, iRow); - const auto& spline = correction.getSpline(iSector, iRow); + auto& info = correction.getRowInfo(iRow); + const auto& spline = correction.getSplineForRow(iRow); auto addVoxel = [&](int iy, int iz, double weight) { auto& vox = vRowVoxels[iy * nZ2Xbins + iz]; @@ -928,6 +902,11 @@ void TPCFastSpaceChargeCorrectionHelper::initInverse(std::vector splineParameters; for (int row = iThread; row < mGeo.getNumberOfRows(); row += mNthreads) { - auto& sectorRowInfo = correction.getSectorRowInfo(sector, row); - sectorRowInfo.gridReal = sectorRowInfo.gridMeasured; + auto& rowInfo = correction.getRowInfo(row); - TPCFastSpaceChargeCorrection::SplineType spline = correction.getSpline(sector, row); + TPCFastSpaceChargeCorrection::SplineType spline = correction.getSplineForRow(row); helper.setSpline(spline, 10, 10); std::vector gridU; @@ -1050,23 +1028,19 @@ void TPCFastSpaceChargeCorrectionHelper::mergeCorrections( auto myThread = [&](int iThread) { for (int row = iThread; row < geo.getNumberOfRows(); row += mNthreads) { - const auto& spline = mainCorrection.getSpline(sector, row); + auto& rowInfo = mainCorrection.getRowInfo(row); + const auto& spline = mainCorrection.getSplineForRow(row); float* splineParameters = mainCorrection.getCorrectionData(sector, row); float* splineParametersInvX = mainCorrection.getCorrectionDataInvX(sector, row); float* splineParametersInvYZ = mainCorrection.getCorrectionDataInvYZ(sector, row); - auto& secRowInfo = mainCorrection.getSectorRowInfo(sector, row); - constexpr int nKnotPar1d = 4; constexpr int nKnotPar2d = nKnotPar1d * 2; constexpr int nKnotPar3d = nKnotPar1d * 3; { // scale the main correction - for (int i = 0; i < 3; i++) { - secRowInfo.maxCorr[i] *= mainScale; - secRowInfo.minCorr[i] *= mainScale; - } + double parscale[4] = {mainScale, mainScale, mainScale, mainScale * mainScale}; for (int iknot = 0, ind = 0; iknot < spline.getNumberOfKnots(); iknot++) { for (int ipar = 0; ipar < nKnotPar1d; ++ipar) { @@ -1099,14 +1073,12 @@ void TPCFastSpaceChargeCorrectionHelper::mergeCorrections( for (int icorr = 0; icorr < additionalCorrections.size(); ++icorr) { const auto& corr = *(additionalCorrections[icorr].first); double scale = additionalCorrections[icorr].second; - auto& linfo = corr.getSectorRowInfo(sector, row); - secRowInfo.updateMaxValues(linfo.getMaxValues(), scale); - secRowInfo.updateMaxValues(linfo.getMinValues(), scale); + auto& linfo = corr.getRowInfo(row); - double scaleU = secRowInfo.gridMeasured.getYscale() / linfo.gridMeasured.getYscale(); - double scaleV = secRowInfo.gridMeasured.getZscale() / linfo.gridMeasured.getZscale(); - double scaleRealU = secRowInfo.gridReal.getYscale() / linfo.gridReal.getYscale(); - double scaleRealV = secRowInfo.gridReal.getZscale() / linfo.gridReal.getZscale(); + double scaleU = rowInfo.gridMeasured.getYscale() / linfo.gridMeasured.getYscale(); + double scaleV = rowInfo.gridMeasured.getZscale() / linfo.gridMeasured.getZscale(); + double scaleRealU = rowInfo.gridReal.getYscale() / linfo.gridReal.getYscale(); + double scaleRealV = rowInfo.gridReal.getZscale() / linfo.gridReal.getZscale(); for (int iu = 0; iu < gridU.getNumberOfKnots(); iu++) { double u = gridU.getKnot(iu).u; @@ -1123,7 +1095,7 @@ void TPCFastSpaceChargeCorrectionHelper::mergeCorrections( corr.convLocalToGrid(sector, row, y, z, lu, lv, ls); ls *= scale; double parscale[4] = {ls, ls * scaleU, ls * scaleV, ls * ls * scaleU * scaleV}; - const auto& spl = corr.getSpline(sector, row); + const auto& spl = corr.getSplineForRow(row); spl.interpolateParametersAtU(corr.getCorrectionData(sector, row), lu, lv, P); for (int ipar = 0, ind = 0; ipar < nKnotPar1d; ++ipar) { for (int idim = 0; idim < 3; idim++, ind++) { @@ -1141,7 +1113,7 @@ void TPCFastSpaceChargeCorrectionHelper::mergeCorrections( double parscale[4] = {ls, ls * scaleRealU, ls * scaleRealV, ls * ls * scaleRealU * scaleRealV}; { // inverse X correction - corr.getSplineInvX(sector, row).interpolateParametersAtU(corr.getCorrectionDataInvX(sector, row), lu, lv, P); + corr.getSplineInvXforRow(row).interpolateParametersAtU(corr.getCorrectionDataInvX(sector, row), lu, lv, P); for (int ipar = 0, ind = 0; ipar < nKnotPar1d; ++ipar) { for (int idim = 0; idim < 1; idim++, ind++) { splineParametersInvX[knotIndex * nKnotPar1d + ind] += parscale[ipar] * P[ind]; @@ -1150,7 +1122,7 @@ void TPCFastSpaceChargeCorrectionHelper::mergeCorrections( } { // inverse YZ correction - corr.getSplineInvYZ(sector, row).interpolateParametersAtU(corr.getCorrectionDataInvYZ(sector, row), lu, lv, P); + corr.getSplineInvYZforRow(row).interpolateParametersAtU(corr.getCorrectionDataInvYZ(sector, row), lu, lv, P); for (int ipar = 0, ind = 0; ipar < nKnotPar1d; ++ipar) { for (int idim = 0; idim < 2; idim++, ind++) { splineParametersInvYZ[knotIndex * nKnotPar2d + ind] += parscale[ipar] * P[ind]; diff --git a/Detectors/TPC/reconstruction/src/TPCFastTransformHelperO2.cxx b/Detectors/TPC/reconstruction/src/TPCFastTransformHelperO2.cxx index 8ae5cbc5fae4a..c08444e0baab3 100644 --- a/Detectors/TPC/reconstruction/src/TPCFastTransformHelperO2.cxx +++ b/Detectors/TPC/reconstruction/src/TPCFastTransformHelperO2.cxx @@ -96,7 +96,7 @@ std::unique_ptr TPCFastTransformHelperO2::create(int64_t TimeS const float t0 = 0.; const float vDrift = 0.f; const long int initTimeStamp = -1; - fastTransform.setCalibration1(initTimeStamp, t0, vDrift); + fastTransform.setCalibration(initTimeStamp, t0, vDrift); fastTransform.finishConstruction(); } @@ -153,7 +153,7 @@ int TPCFastTransformHelperO2::updateCalibrationImpl(T& fastTransform, int64_t Ti const double t0 = (driftTimeOffset + elParam.getAverageShapingTime()) / elParam.ZbinWidth; - fastTransform.setCalibration1(TimeStamp, t0, vDrift); + fastTransform.setCalibration(TimeStamp, t0, vDrift); return 0; } diff --git a/GPU/TPCFastTransformation/CorrectionMapsHelper.h b/GPU/TPCFastTransformation/CorrectionMapsHelper.h index 39c5ffc73b1da..095bb837eaacd 100644 --- a/GPU/TPCFastTransformation/CorrectionMapsHelper.h +++ b/GPU/TPCFastTransformation/CorrectionMapsHelper.h @@ -122,7 +122,7 @@ class CorrectionMapsHelper { if (mCorrMapMShape) { // just check for the first spline the number of knots which are 4 in case of default spline object - return mCorrMapMShape->getCorrection().getSpline(0, 0).getNumberOfKnots() == 4; + return mCorrMapMShape->getCorrection().getSplineForRow(0).getNumberOfKnots() == 4; } return true; } diff --git a/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.cxx b/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.cxx index 49305e3ed7909..7632c3587ecc8 100644 --- a/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.cxx +++ b/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.cxx @@ -34,7 +34,7 @@ TPCFastSpaceChargeCorrection::TPCFastSpaceChargeCorrection() mScenarioPtr(nullptr), mTimeStamp(-1), mCorrectionData{nullptr, nullptr, nullptr}, - mCorrectionDataSize{0, 0, 0} + mSectorDataSizeBytes{0, 0, 0} { // Default Constructor: creates an empty uninitialized object } @@ -61,7 +61,7 @@ void TPCFastSpaceChargeCorrection::destroy() mTimeStamp = -1; for (int32_t is = 0; is < 3; is++) { mCorrectionData[is] = nullptr; - mCorrectionDataSize[is] = 0; + mSectorDataSizeBytes[is] = 0; } FlatObject::destroy(); } @@ -98,9 +98,9 @@ void TPCFastSpaceChargeCorrection::cloneFromObject(const TPCFastSpaceChargeCorre mTimeStamp = obj.mTimeStamp; - mCorrectionDataSize[0] = obj.mCorrectionDataSize[0]; - mCorrectionDataSize[1] = obj.mCorrectionDataSize[1]; - mCorrectionDataSize[2] = obj.mCorrectionDataSize[2]; + mSectorDataSizeBytes[0] = obj.mSectorDataSizeBytes[0]; + mSectorDataSizeBytes[1] = obj.mSectorDataSizeBytes[1]; + mSectorDataSizeBytes[2] = obj.mSectorDataSizeBytes[2]; // variable-size data mScenarioPtr = obj.mScenarioPtr; @@ -110,8 +110,8 @@ void TPCFastSpaceChargeCorrection::cloneFromObject(const TPCFastSpaceChargeCorre mClassVersion = obj.mClassVersion; - for (int32_t i = 0; i < TPCFastTransformGeo::getNumberOfSectors() * TPCFastTransformGeo::getMaxNumberOfRows(); i++) { - mSectorRowInfos[i] = obj.mSectorRowInfos[i]; + for (int32_t i = 0; i < TPCFastTransformGeo::getMaxNumberOfRows(); i++) { + mRowInfos[i] = obj.mRowInfos[i]; } relocateBufferPointers(oldFlatBufferPtr, mFlatBufferPtr); @@ -126,37 +126,6 @@ void TPCFastSpaceChargeCorrection::moveBufferTo(char* newFlatBufferPtr) relocateBufferPointers(oldFlatBufferPtr, mFlatBufferPtr); } -void TPCFastSpaceChargeCorrection::setActualBufferAddressOld(char* actualFlatBufferPtr) -{ - /// Sets the actual location of the external flat buffer after it has been moved (e.g. to another maschine) - - if (mClassVersion != 4) { - LOG(error) << "TPCFastSpaceChargeCorrection::setActualBufferAddress() called with class version " << mClassVersion << ". This is not supported."; - return; - } - - FlatObject::setActualBufferAddress(actualFlatBufferPtr); - - size_t scSize = sizeof(SplineType) * mNumberOfScenarios; - - mScenarioPtr = reinterpret_cast(mFlatBufferPtr); - - size_t scBufferOffset = alignSize(scSize, SplineType::getBufferAlignmentBytes()); - size_t scBufferSize = 0; - - for (int32_t i = 0; i < mNumberOfScenarios; i++) { - SplineType& sp = mScenarioPtr[i]; - sp.setActualBufferAddress(mFlatBufferPtr + scBufferOffset + scBufferSize); - scBufferSize = alignSize(scBufferSize + sp.getFlatBufferSize(), sp.getBufferAlignmentBytes()); - } - size_t bufferSize = scBufferOffset + scBufferSize; - for (int32_t is = 0; is < 3; is++) { - size_t correctionDataOffset = alignSize(bufferSize, SplineType::getParameterAlignmentBytes()); - mCorrectionData[is] = reinterpret_cast(mFlatBufferPtr + correctionDataOffset); - bufferSize = correctionDataOffset + mCorrectionDataSize[is]; - } -} - void TPCFastSpaceChargeCorrection::setActualBufferAddress(char* actualFlatBufferPtr) { /// Sets the actual location of the external flat buffer after it has been moved (e.g. to another maschine) @@ -178,9 +147,8 @@ void TPCFastSpaceChargeCorrection::setActualBufferAddress(char* actualFlatBuffer } size_t bufferSize = scBufferOffset + scBufferSize; for (int32_t is = 0; is < 3; is++) { - size_t correctionDataOffset = alignSize(bufferSize, SplineType::getParameterAlignmentBytes()); - mCorrectionData[is] = reinterpret_cast(mFlatBufferPtr + correctionDataOffset); - bufferSize = correctionDataOffset + mCorrectionDataSize[is]; + mCorrectionData[is] = reinterpret_cast(mFlatBufferPtr + bufferSize); + bufferSize += mSectorDataSizeBytes[is] * mGeo.getNumberOfSectors(); } return; } @@ -232,59 +200,31 @@ void TPCFastSpaceChargeCorrection::setActualBufferAddress(char* actualFlatBuffer auto* oldRowInfos = reinterpret_cast(mFlatBufferPtr + oldRowsOffset); auto* oldSectorRowInfos = reinterpret_cast(mFlatBufferPtr + oldSectorRowsOffset); - size_t sectorDataSize[3]; - for (int32_t is = 0; is < 3; is++) { - sectorDataSize[is] = mCorrectionDataSize[is] / mGeo.getNumberOfSectors(); - } + int32_t iSector = 0; - for (int32_t iSector = 0; iSector < mGeo.getNumberOfSectors(); iSector++) { - - for (int32_t iRow = 0; iRow < mGeo.getNumberOfRows(); iRow++) { - RowInfoVersion3& oldRowInfo = oldRowInfos[iRow]; - SectorRowInfoVersion3& oldSectorRowInfo = oldSectorRowInfos[mGeo.getNumberOfRows() * iSector + iRow]; - - // the spline buffer is not yet initialised, don't try to access knot positions etc - const auto& spline = oldScenarioPtr[oldRowInfo.splineScenarioID]; - - SectorRowInfo& newSectorRow = getSectorRowInfo(iSector, iRow); + for (int32_t iRow = 0; iRow < mGeo.getNumberOfRows(); iRow++) { + RowInfoVersion3& oldRowInfo = oldRowInfos[iRow]; + SectorRowInfoVersion3& oldSectorRowInfo = oldSectorRowInfos[mGeo.getNumberOfRows() * iSector + iRow]; - newSectorRow.splineScenarioID = oldRowInfo.splineScenarioID; - for (int32_t is = 0; is < 3; is++) { - newSectorRow.dataOffsetBytes[is] = sectorDataSize[is] * iSector + oldRowInfo.dataOffsetBytes[is]; - } + // the spline buffer is not yet initialised, don't try to access knot positions etc + const auto& spline = oldScenarioPtr[oldRowInfo.splineScenarioID]; - { // grid for the measured coordinates - float y0 = mGeo.getRowInfo(iRow).yMin; - float yScale = spline.getGridX1().getUmax() / mGeo.getRowInfo(iRow).getYwidth(); - float zReadout = mGeo.getZreadout(iSector); - float zOut = mGeo.getTPCzLength() - oldSectorRowInfo.gridV0; - float z0 = -3.; - float zScale = spline.getGridX2().getUmax() / (zOut - z0); - if (iSector >= mGeo.getNumberOfSectorsA()) { - zOut = -zOut; - z0 = zOut; - } - newSectorRow.gridMeasured.set(y0, yScale, z0, zScale, zOut, zReadout); - } + RowInfo& newRowInfo = getRowInfo(iRow); - { // grid for the real coordinates - float y0 = oldSectorRowInfo.gridCorrU0; - float yScale = oldSectorRowInfo.scaleCorrUtoGrid; - float zReadout = mGeo.getZreadout(iSector); - float zOut = mGeo.getTPCzLength() - oldSectorRowInfo.gridCorrV0; - float zScale = oldSectorRowInfo.scaleCorrVtoGrid; - float z0 = zOut - spline.getGridX2().getUmax() / zScale; - if (iSector >= mGeo.getNumberOfSectorsA()) { - zOut = -zOut; - z0 = zOut; - } - newSectorRow.gridReal.set(y0, yScale, z0, zScale, zOut, zReadout); - } + newRowInfo.splineScenarioID = oldRowInfo.splineScenarioID; + for (int32_t is = 0; is < 3; is++) { + newRowInfo.dataOffsetBytes[is] = oldRowInfo.dataOffsetBytes[is]; + } - newSectorRow.resetMaxValues(); - newSectorRow.updateMaxValues(-100.f, -100.f, -100.f); - newSectorRow.updateMaxValues(100.f, 100.f, 100.f); + { // grid for the measured coordinates + float y0 = mGeo.getRowInfo(iRow).yMin; + float yScale = spline.getGridX1().getUmax() / mGeo.getRowInfo(iRow).getYwidth(); + float zOut = mGeo.getTPCzLength() - oldSectorRowInfo.gridV0; + float z0 = -3.; + float zScale = spline.getGridX2().getUmax() / (zOut - z0); + newRowInfo.gridMeasured.set(y0, yScale, z0, zScale, zOut, mGeo.getTPCzLength()); } + newRowInfo.gridReal = newRowInfo.gridMeasured; } } @@ -315,11 +255,11 @@ void TPCFastSpaceChargeCorrection::setActualBufferAddress(char* actualFlatBuffer for (int32_t is = 0; is < 3; is++) { size_t oldCorrectionDataOffset = alignSize(oldBufferSize, SplineType::getParameterAlignmentBytes()); - size_t correctionDataOffset = alignSize(bufferSize, SplineType::getParameterAlignmentBytes()); + size_t correctionDataOffset = bufferSize; mCorrectionData[is] = reinterpret_cast(mFlatBufferPtr + correctionDataOffset); - memmove(mCorrectionData[is], mFlatBufferPtr + oldCorrectionDataOffset, mCorrectionDataSize[is]); - oldBufferSize = oldCorrectionDataOffset + mCorrectionDataSize[is]; - bufferSize = correctionDataOffset + mCorrectionDataSize[is]; + memmove(mCorrectionData[is], mFlatBufferPtr + oldCorrectionDataOffset, mSectorDataSizeBytes[is] * mGeo.getNumberOfSectors()); + oldBufferSize = oldCorrectionDataOffset + mSectorDataSizeBytes[is] * mGeo.getNumberOfSectors(); + bufferSize = correctionDataOffset + mSectorDataSizeBytes[is] * mGeo.getNumberOfSectors(); } mFlatBufferSize = bufferSize; @@ -329,8 +269,8 @@ void TPCFastSpaceChargeCorrection::setActualBufferAddress(char* actualFlatBuffer bool isAside = (iSector < mGeo.getNumberOfSectorsA()); for (int32_t iRow = 0; iRow < mGeo.getNumberOfRows(); iRow++) { - SectorRowInfo& sectorRow = getSectorRowInfo(iSector, iRow); - const auto& spline = mScenarioPtr[sectorRow.splineScenarioID]; + RowInfo& rowInfo = getRowInfo(iRow); + const auto& spline = mScenarioPtr[rowInfo.splineScenarioID]; int nSplineDimensions[3] = {3, 1, 2}; @@ -348,17 +288,17 @@ void TPCFastSpaceChargeCorrection::setActualBufferAddress(char* actualFlatBuffer } }; - // reorder knots for the A side Y == old U, Z == - old V + // reorder knots for the A side U == old U, V == - old V if (isAside) { for (int32_t i = 0; i < spline.getGridX1().getNumberOfKnots(); i++) { for (int32_t j = 0; j < spline.getGridX2().getNumberOfKnots() / 2; j++) { swapKnots(i, j, i, spline.getGridX2().getNumberOfKnots() - 1 - j); } } - } else { // reorder knots for the C side Y == - old U, Z == old V + } else { // reorder knots for the C side U == - old U, V == - old V for (int32_t i = 0; i < spline.getGridX1().getNumberOfKnots() / 2; i++) { for (int32_t j = 0; j < spline.getGridX2().getNumberOfKnots(); j++) { - swapKnots(i, j, spline.getGridX1().getNumberOfKnots() - 1 - i, j); + swapKnots(i, j, spline.getGridX1().getNumberOfKnots() - 1 - i, spline.getGridX2().getNumberOfKnots() - 1 - j); } } } @@ -370,10 +310,11 @@ void TPCFastSpaceChargeCorrection::setActualBufferAddress(char* actualFlatBuffer for (int iDim = 0; iDim < nDim; iDim++) { if (isAside) { data[nKnotParameters * iKnot + nDim * 1 + iDim] *= -1; // invert Z derivatives on A side + data[nKnotParameters * iKnot + nDim * 3 + iDim] *= -1; // invert cross derivatives on A side } else { + data[nKnotParameters * iKnot + nDim * 1 + iDim] *= -1; // invert Z derivatives on C side data[nKnotParameters * iKnot + nDim * 2 + iDim] *= -1; // invert Y derivatives on C side } - data[nKnotParameters * iKnot + nDim * 3 + iDim] *= -1; // invert cross derivatives on both sides } // new correction directions if (iSpline == 0) { // dX,dU,dV -> dX,dY,dZ @@ -442,7 +383,7 @@ void TPCFastSpaceChargeCorrection::print() const mGeo.print(); LOG(info) << " mNumberOfScenarios = " << mNumberOfScenarios; LOG(info) << " mTimeStamp = " << mTimeStamp; - LOG(info) << " mCorrectionDataSize = " << mCorrectionDataSize[0] << " " << mCorrectionDataSize[1] << " " << mCorrectionDataSize[2]; + LOG(info) << " mSectorDataSizeBytes = " << mSectorDataSizeBytes[0] << " " << mSectorDataSizeBytes[1] << " " << mSectorDataSizeBytes[2]; if (mScenarioPtr) { for (int32_t i = 0; i < mNumberOfScenarios; i++) { @@ -455,7 +396,7 @@ void TPCFastSpaceChargeCorrection::print() const for (int32_t is = 0; is < mGeo.getNumberOfSectors(); is++) { for (int32_t ir = 0; ir < mGeo.getNumberOfRows(); ir++) { LOG(info) << "sector " << is << " row " << ir << ": "; - const SplineType& spline = getSpline(is, ir); + const SplineType& spline = getSplineForRow(ir); const float* d = getCorrectionData(is, ir); int32_t k = 0; for (int32_t i = 0; i < spline.getGridX1().getNumberOfKnots(); i++) { @@ -488,22 +429,14 @@ void TPCFastSpaceChargeCorrection::startConstruction(const TPCFastTransformGeo& assert(mConstructionScenarios != nullptr); - for (int32_t i = 0; i < mGeo.getNumberOfSectors(); i++) { - for (int32_t j = 0; j < mGeo.getNumberOfRows(); j++) { - auto& row = mSectorRowInfos[mGeo.getMaxNumberOfRows() * i + j]; - row.splineScenarioID = -1; - row.gridReal = {}; - row.gridMeasured = {}; - row.dataOffsetBytes[0] = 0; - row.dataOffsetBytes[1] = 0; - row.dataOffsetBytes[2] = 0; - row.minCorr[0] = 0; - row.minCorr[1] = 0; - row.minCorr[2] = 0; - row.maxCorr[0] = 0; - row.maxCorr[1] = 0; - row.maxCorr[2] = 0; - } + for (int32_t i = 0; i < mGeo.getNumberOfRows(); i++) { + auto& row = mRowInfos[i]; + row.splineScenarioID = -1; + row.gridReal = {}; + row.gridMeasured = {}; + row.dataOffsetBytes[0] = 0; + row.dataOffsetBytes[1] = 0; + row.dataOffsetBytes[2] = 0; } for (int32_t i = 0; i < mNumberOfScenarios; i++) { @@ -515,18 +448,18 @@ void TPCFastSpaceChargeCorrection::startConstruction(const TPCFastTransformGeo& mScenarioPtr = nullptr; for (int32_t s = 0; s < 3; s++) { mCorrectionData[s] = nullptr; - mCorrectionDataSize[s] = 0; + mSectorDataSizeBytes[s] = 0; } mClassVersion = 4; } -void TPCFastSpaceChargeCorrection::setRowScenarioID(int32_t iSector, int32_t iRow, int32_t iScenario) +void TPCFastSpaceChargeCorrection::setRowScenarioID(int32_t iRow, int32_t iScenario) { /// Initializes a TPC row assert(mConstructionMask & ConstructionState::InProgress); assert(iSector >= 0 && iSector < mGeo.getNumberOfSectors()); assert(iRow >= 0 && iRow < mGeo.getNumberOfRows() && iScenario >= 0 && iScenario < mNumberOfScenarios); - auto& row = getSectorRowInfo(iSector, iRow); + auto& row = getRowInfo(iRow); row.splineScenarioID = iScenario; for (int32_t s = 0; s < 3; s++) { row.dataOffsetBytes[s] = 0; @@ -549,12 +482,10 @@ void TPCFastSpaceChargeCorrection::finishConstruction() assert(mConstructionMask & ConstructionState::InProgress); - for (int32_t i = 0; i < mGeo.getNumberOfSectors(); i++) { - for (int32_t j = 0; j < mGeo.getNumberOfRows(); j++) { - [[maybe_unused]] SectorRowInfo& row = getSectorRowInfo(i, j); - assert(row.splineScenarioID >= 0); - assert(row.splineScenarioID < mNumberOfScenarios); - } + for (int32_t j = 0; j < mGeo.getNumberOfRows(); j++) { + [[maybe_unused]] RowInfo& row = getRowInfo(j); + assert(row.splineScenarioID >= 0); + assert(row.splineScenarioID < mNumberOfScenarios); } for (int32_t i = 0; i < mNumberOfScenarios; i++) { @@ -578,18 +509,15 @@ void TPCFastSpaceChargeCorrection::finishConstruction() size_t bufferSize = scBufferOffsets[0] + scBufferSize; size_t correctionDataOffset[3]; for (int32_t is = 0; is < 3; is++) { - correctionDataOffset[is] = alignSize(bufferSize, SplineType::getParameterAlignmentBytes()); - mCorrectionDataSize[is] = 0; - for (int32_t i = 0; i < mGeo.getNumberOfSectors(); i++) { - for (int32_t j = 0; j < mGeo.getNumberOfRows(); j++) { - SectorRowInfo& row = getSectorRowInfo(i, j); - SplineType& spline = mConstructionScenarios[row.splineScenarioID]; - row.dataOffsetBytes[is] = alignSize(mCorrectionDataSize[is], SplineType::getParameterAlignmentBytes()); - mCorrectionDataSize[is] = row.dataOffsetBytes[is] + spline.getSizeOfParameters(); - } + correctionDataOffset[is] = bufferSize; + mSectorDataSizeBytes[is] = 0; + for (int32_t j = 0; j < mGeo.getNumberOfRows(); j++) { + RowInfo& row = getRowInfo(j); + SplineType& spline = mConstructionScenarios[row.splineScenarioID]; + row.dataOffsetBytes[is] = mSectorDataSizeBytes[is]; + mSectorDataSizeBytes[is] += spline.getSizeOfParameters(); } - mCorrectionDataSize[is] = alignSize(mCorrectionDataSize[is], SplineType::getParameterAlignmentBytes()); - bufferSize = correctionDataOffset[is] + mCorrectionDataSize[is]; + bufferSize += mSectorDataSizeBytes[is] * mGeo.getNumberOfSectors(); } FlatObject::finishConstruction(bufferSize); @@ -616,11 +544,22 @@ void TPCFastSpaceChargeCorrection::finishConstruction() GPUd() void TPCFastSpaceChargeCorrection::setNoCorrection() { // initialise all corrections to 0. + + for (int32_t row = 0; row < mGeo.getNumberOfRows(); row++) { + const SplineType& spline = getSplineForRow(row); + RowInfo& info = getRowInfo(row); + float y0 = mGeo.getRowInfo(row).getYmin(); + float yScale = spline.getGridX1().getUmax() / mGeo.getRowInfo(row).getYwidth(); + float z0 = 0.; + float zScale = spline.getGridX2().getUmax() / mGeo.getTPCzLength(); + info.gridMeasured.set(y0, yScale, z0, zScale, mGeo.getTPCzLength(), mGeo.getTPCzLength()); + info.gridReal = info.gridMeasured; + } // row + for (int32_t sector = 0; sector < mGeo.getNumberOfSectors(); sector++) { for (int32_t row = 0; row < mGeo.getNumberOfRows(); row++) { - const SplineType& spline = getSpline(sector, row); - + const SplineType& spline = getSplineForRow(row); for (int32_t is = 0; is < 3; is++) { float* data = getCorrectionData(sector, row, is); int32_t nPar = spline.getNumberOfParameters(); @@ -634,17 +573,6 @@ GPUd() void TPCFastSpaceChargeCorrection::setNoCorrection() data[i] = 0.f; } } - - SectorRowInfo& info = getSectorRowInfo(sector, row); - - float y0 = mGeo.getRowInfo(row).getYmin(); - float yScale = spline.getGridX1().getUmax() / mGeo.getRowInfo(row).getYwidth(); - float z0 = mGeo.getZmin(sector); - float zScale = spline.getGridX2().getUmax() / mGeo.getTPCzLength(); - float zReadout = mGeo.getZreadout(sector); - info.gridMeasured.set(y0, yScale, z0, zScale, zReadout, zReadout); - - info.gridReal = info.gridMeasured; } // row } // sector } @@ -653,10 +581,8 @@ void TPCFastSpaceChargeCorrection::constructWithNoCorrection(const TPCFastTransf { const int32_t nCorrectionScenarios = 1; startConstruction(geo, nCorrectionScenarios); - for (int32_t sector = 0; sector < geo.getNumberOfSectors(); sector++) { - for (int32_t row = 0; row < geo.getNumberOfRows(); row++) { - setRowScenarioID(sector, row, 0); - } + for (int32_t row = 0; row < geo.getNumberOfRows(); row++) { + setRowScenarioID(row, 0); } { TPCFastSpaceChargeCorrection::SplineType spline; diff --git a/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h b/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h index a32c835ad7731..09704bb5706e1 100644 --- a/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h +++ b/GPU/TPCFastTransformation/TPCFastSpaceChargeCorrection.h @@ -42,7 +42,7 @@ class TPCFastSpaceChargeCorrection : public FlatObject friend class TPCFastTransformPOD; public: - // obsolete structure, declared here only for backward compatibility + // obsolete structure, declared here only for the backward compatibility struct SliceInfo { ClassDefNV(SliceInfo, 2); }; @@ -79,75 +79,36 @@ class TPCFastSpaceChargeCorrection : public FlatObject } /// convert local y, z to internal grid coordinates u,v, and spline scale - GPUdi() void convLocalToGridUntruncated(float y, float z, float& u, float& v, float& s) const + GPUdi() void convLocalToGridUntruncated(int sector, float y, float z, float& u, float& v, float& s) const { + if (sector >= TPCFastTransformGeo::getNumberOfSectorsA()) { + z = -z; + } u = (y - y0) * yScale; v = (z - z0) * zScale; s = getSpineScaleForZ(z); } /// convert internal grid coordinates u,v to local y, z - GPUdi() void convGridToLocal(float gridU, float gridV, float& y, float& z) const + GPUdi() void convGridToLocal(int sector, float gridU, float gridV, float& y, float& z) const { y = y0 + gridU / yScale; z = z0 + gridV / zScale; + if (sector >= TPCFastTransformGeo::getNumberOfSectorsA()) { + z = -z; + } } ClassDefNV(GridInfo, 1); }; - struct SectorRowInfo { + struct RowInfo { int32_t splineScenarioID{0}; ///< scenario index (which of Spline2D splines to use) size_t dataOffsetBytes[3]{0}; ///< offset for the spline data withing a TPC sector - GridInfo gridMeasured; ///< grid info for measured coordinates - GridInfo gridReal; ///< grid info for real coordinates + GridInfo gridMeasured; ///< grid info for the measured coordinates + GridInfo gridReal; ///< grid info for the real coordinates - float minCorr[3]{-10.f, -10.f, -10.f}; ///< min correction for dX, dY, dZ - float maxCorr[3]{10.f, 10.f, 10.f}; ///< max correction for dX, dY, dZ - - void resetMaxValues() - { - minCorr[0] = -1.f; - maxCorr[0] = 1.f; - minCorr[1] = -1.f; - maxCorr[1] = 1.f; - minCorr[2] = -1.f; - maxCorr[2] = 1.f; - } - - void updateMaxValues(float dx, float du, float dv) - { - minCorr[0] = GPUCommonMath::Min(minCorr[0], dx); - maxCorr[0] = GPUCommonMath::Max(maxCorr[0], dx); - - minCorr[1] = GPUCommonMath::Min(minCorr[1], du); - maxCorr[1] = GPUCommonMath::Max(maxCorr[1], du); - - minCorr[2] = GPUCommonMath::Min(minCorr[2], dv); - maxCorr[2] = GPUCommonMath::Max(maxCorr[2], dv); - } - -#ifndef GPUCA_GPUCODE_DEVICE - void updateMaxValues(std::array dxdudv, float scale) - { - float dx = dxdudv[0] * scale; - float du = dxdudv[1] * scale; - float dv = dxdudv[2] * scale; - updateMaxValues(dx, du, dv); - } - - std::array getMaxValues() const - { - return {maxCorr[0], maxCorr[1], maxCorr[2]}; - } - - std::array getMinValues() const - { - return {minCorr[0], minCorr[1], minCorr[2]}; - } -#endif - - ClassDefNV(SectorRowInfo, 2); + ClassDefNV(RowInfo, 2); }; typedef Spline2D SplineTypeXYZ; @@ -195,7 +156,6 @@ class TPCFastSpaceChargeCorrection : public FlatObject /// Moving the class with its external buffer to another location - void setActualBufferAddressOld(char* actualFlatBufferPtr); void setActualBufferAddress(char* actualFlatBufferPtr); void setFutureBufferAddress(char* futureFlatBufferPtr); @@ -205,7 +165,7 @@ class TPCFastSpaceChargeCorrection : public FlatObject void startConstruction(const TPCFastTransformGeo& geo, int32_t numberOfSplineScenarios); /// Initializes a TPC row - void setRowScenarioID(int32_t iSector, int32_t iRow, int32_t iScenario); + void setRowScenarioID(int32_t iRow, int32_t iScenario); /// Sets approximation scenario void setSplineScenario(int32_t scenarioIndex, const SplineType& spline); @@ -224,10 +184,10 @@ class TPCFastSpaceChargeCorrection : public FlatObject GPUdi() void setTimeStamp(int64_t v) { mTimeStamp = v; } /// Gives const pointer to a spline - GPUd() const SplineType& getSpline(int32_t sector, int32_t row) const; + GPUd() const SplineType& getSplineForRow(int32_t row) const; /// Gives pointer to a spline - GPUd() SplineType& getSpline(int32_t sector, int32_t row); + GPUd() SplineType& getSplineForRow(int32_t row); /// Gives pointer to spline data GPUd() float* getCorrectionData(int32_t sector, int32_t row, int32_t iSpline = 0); @@ -236,10 +196,10 @@ class TPCFastSpaceChargeCorrection : public FlatObject GPUd() const float* getCorrectionData(int32_t sector, int32_t row, int32_t iSpline = 0) const; /// Gives const pointer to a spline for the inverse X correction - GPUd() const SplineTypeInvX& getSplineInvX(int32_t sector, int32_t row) const; + GPUd() const SplineTypeInvX& getSplineInvXforRow(int32_t row) const; /// Gives pointer to a spline for the inverse X correction - GPUd() SplineTypeInvX& getSplineInvX(int32_t sector, int32_t row); + GPUd() SplineTypeInvX& getSplineInvXforRow(int32_t row); /// Gives pointer to spline data for the inverse X correction GPUd() float* getCorrectionDataInvX(int32_t sector, int32_t row); @@ -248,10 +208,10 @@ class TPCFastSpaceChargeCorrection : public FlatObject GPUd() const float* getCorrectionDataInvX(int32_t sector, int32_t row) const; /// Gives const pointer to a spline for the inverse YZ correction - GPUd() const SplineTypeInvYZ& getSplineInvYZ(int32_t sector, int32_t row) const; + GPUd() const SplineTypeInvYZ& getSplineInvYZforRow(int32_t row) const; /// Gives pointer to a spline for the inverse YZ correction - GPUd() SplineTypeInvYZ& getSplineInvYZ(int32_t sector, int32_t row); + GPUd() SplineTypeInvYZ& getSplineInvYZforRow(int32_t row); /// Gives pointer to spline data for the inverse YZ correction GPUd() float* getCorrectionDataInvYZ(int32_t sector, int32_t row); @@ -301,16 +261,16 @@ class TPCFastSpaceChargeCorrection : public FlatObject /// Gives the time stamp of the current calibaration parameters int64_t getTimeStamp() const { return mTimeStamp; } - /// Gives TPC sector & row info - GPUdi() const SectorRowInfo& getSectorRowInfo(int32_t sector, int32_t row) const + /// Gives TPC row info + GPUdi() const RowInfo& getRowInfo(int32_t row) const { - return mSectorRowInfos[mGeo.getMaxNumberOfRows() * sector + row]; + return mRowInfos[row]; } - /// Gives TPC sector & row info - GPUdi() SectorRowInfo& getSectorRowInfo(int32_t sector, int32_t row) + /// Gives TPC row info + GPUdi() RowInfo& getRowInfo(int32_t row) { - return mSectorRowInfos[mGeo.getMaxNumberOfRows() * sector + row]; + return mRowInfos[row]; } #if !defined(GPUCA_GPUCODE) @@ -325,6 +285,8 @@ class TPCFastSpaceChargeCorrection : public FlatObject /// release temporary memory used during construction void releaseConstructionMemory(); + static constexpr float kMaxCorrection = 100.f; ///< maximum correction value, used to protect from FPEs + /// _______________ Data members _______________________________________________ /// _______________ Construction control _______________________________________________ @@ -345,14 +307,14 @@ class TPCFastSpaceChargeCorrection : public FlatObject char* mCorrectionData[3]; //! (transient!!) pointer to the spline data in the flat buffer - size_t mCorrectionDataSize[3]; ///< size of the data per transformation (direct, inverseX, inverse YZ) in the flat buffer + size_t mSectorDataSizeBytes[3]; ///< size of the sector data per transformation (direct, inverseX, inverse YZ) in the flat buffer /// Class version. It is used to read older versions from disc. /// The default version 3 is the one before this field was introduced. /// The actual version must be set in startConstruction(). int32_t mClassVersion{3}; - SectorRowInfo mSectorRowInfos[TPCFastTransformGeo::getNumberOfSectors() * TPCFastTransformGeo::getMaxNumberOfRows()]; ///< SectorRowInfo array + RowInfo mRowInfos[TPCFastTransformGeo::getMaxNumberOfRows()]; ///< RowInfo array ClassDefNV(TPCFastSpaceChargeCorrection, 4); }; @@ -361,40 +323,42 @@ class TPCFastSpaceChargeCorrection : public FlatObject /// Inline implementations of some methods /// ==================================================== -GPUdi() const TPCFastSpaceChargeCorrection::SplineType& TPCFastSpaceChargeCorrection::getSpline(int32_t sector, int32_t row) const +GPUdi() const TPCFastSpaceChargeCorrection::SplineType& TPCFastSpaceChargeCorrection::getSplineForRow(int32_t row) const { /// Gives const pointer to spline - return mScenarioPtr[getSectorRowInfo(sector, row).splineScenarioID]; + return mScenarioPtr[getRowInfo(row).splineScenarioID]; } -GPUdi() TPCFastSpaceChargeCorrection::SplineType& TPCFastSpaceChargeCorrection::getSpline(int32_t sector, int32_t row) +GPUdi() TPCFastSpaceChargeCorrection::SplineType& TPCFastSpaceChargeCorrection::getSplineForRow(int32_t row) { /// Gives pointer to spline - return mScenarioPtr[getSectorRowInfo(sector, row).splineScenarioID]; + return mScenarioPtr[getRowInfo(row).splineScenarioID]; } GPUdi() float* TPCFastSpaceChargeCorrection::getCorrectionData(int32_t sector, int32_t row, int32_t iSpline) { /// Gives pointer to spline data - return reinterpret_cast(mCorrectionData[iSpline] + getSectorRowInfo(sector, row).dataOffsetBytes[iSpline]); + size_t offset = sector * mSectorDataSizeBytes[iSpline] + getRowInfo(row).dataOffsetBytes[iSpline]; + return reinterpret_cast(mCorrectionData[iSpline] + offset); } GPUdi() const float* TPCFastSpaceChargeCorrection::getCorrectionData(int32_t sector, int32_t row, int32_t iSpline) const { /// Gives pointer to spline data - return reinterpret_cast(mCorrectionData[iSpline] + getSectorRowInfo(sector, row).dataOffsetBytes[iSpline]); + size_t offset = sector * mSectorDataSizeBytes[iSpline] + getRowInfo(row).dataOffsetBytes[iSpline]; + return reinterpret_cast(mCorrectionData[iSpline] + offset); } -GPUdi() TPCFastSpaceChargeCorrection::SplineTypeInvX& TPCFastSpaceChargeCorrection::getSplineInvX(int32_t sector, int32_t row) +GPUdi() TPCFastSpaceChargeCorrection::SplineTypeInvX& TPCFastSpaceChargeCorrection::getSplineInvXforRow(int32_t row) { /// Gives pointer to spline for the inverse X correction - return reinterpret_cast(getSpline(sector, row)); + return reinterpret_cast(getSplineForRow(row)); } -GPUdi() const TPCFastSpaceChargeCorrection::SplineTypeInvX& TPCFastSpaceChargeCorrection::getSplineInvX(int32_t sector, int32_t row) const +GPUdi() const TPCFastSpaceChargeCorrection::SplineTypeInvX& TPCFastSpaceChargeCorrection::getSplineInvXforRow(int32_t row) const { /// Gives const pointer to spline for the inverse X correction - return reinterpret_cast(getSpline(sector, row)); + return reinterpret_cast(getSplineForRow(row)); } GPUdi() float* TPCFastSpaceChargeCorrection::getCorrectionDataInvX(int32_t sector, int32_t row) @@ -409,16 +373,16 @@ GPUdi() const float* TPCFastSpaceChargeCorrection::getCorrectionDataInvX(int32_t return getCorrectionData(sector, row, 1); } -GPUdi() TPCFastSpaceChargeCorrection::SplineTypeInvYZ& TPCFastSpaceChargeCorrection::getSplineInvYZ(int32_t sector, int32_t row) +GPUdi() TPCFastSpaceChargeCorrection::SplineTypeInvYZ& TPCFastSpaceChargeCorrection::getSplineInvYZforRow(int32_t row) { /// Gives pointer to spline for the inverse YZ correction - return reinterpret_cast(getSpline(sector, row)); + return reinterpret_cast(getSplineForRow(row)); } -GPUdi() const TPCFastSpaceChargeCorrection::SplineTypeInvYZ& TPCFastSpaceChargeCorrection::getSplineInvYZ(int32_t sector, int32_t row) const +GPUdi() const TPCFastSpaceChargeCorrection::SplineTypeInvYZ& TPCFastSpaceChargeCorrection::getSplineInvYZforRow(int32_t row) const { /// Gives const pointer to spline for the inverse YZ correction - return reinterpret_cast(getSpline(sector, row)); + return reinterpret_cast(getSplineForRow(row)); } GPUdi() float* TPCFastSpaceChargeCorrection::getCorrectionDataInvYZ(int32_t sector, int32_t row) @@ -437,8 +401,8 @@ GPUdi() void TPCFastSpaceChargeCorrection::convLocalToGrid(int32_t sector, int32 { /// convert local y, z to internal grid coordinates u,v /// return values: u, v, scaling factor - const SplineType& spline = getSpline(sector, row); - getSectorRowInfo(sector, row).gridMeasured.convLocalToGridUntruncated(y, z, u, v, s); + const SplineType& spline = getSplineForRow(row); + getRowInfo(row).gridMeasured.convLocalToGridUntruncated(sector, y, z, u, v, s); // shrink to the grid u = GPUCommonMath::Clamp(u, 0.f, (float)spline.getGridX1().getUmax()); v = GPUCommonMath::Clamp(v, 0.f, (float)spline.getGridX2().getUmax()); @@ -448,8 +412,8 @@ GPUdi() bool TPCFastSpaceChargeCorrection::isLocalInsideGrid(int32_t sector, int { /// check if local y, z are inside the grid float u, v, s; - getSectorRowInfo(sector, row).gridMeasured.convLocalToGridUntruncated(y, z, u, v, s); - const auto& spline = getSpline(sector, row); + getRowInfo(row).gridMeasured.convLocalToGridUntruncated(sector, y, z, u, v, s); + const auto& spline = getSplineForRow(row); // shrink to the grid if (u < 0.f || u > (float)spline.getGridX1().getUmax() || // v < 0.f || v > (float)spline.getGridX2().getUmax()) { @@ -462,8 +426,8 @@ GPUdi() bool TPCFastSpaceChargeCorrection::isRealLocalInsideGrid(int32_t sector, { /// check if local y, z are inside the grid float u, v, s; - getSectorRowInfo(sector, row).gridReal.convLocalToGridUntruncated(y, z, u, v, s); - const auto& spline = getSpline(sector, row); + getRowInfo(row).gridReal.convLocalToGridUntruncated(sector, y, z, u, v, s); + const auto& spline = getSplineForRow(row); // shrink to the grid if (u < 0.f || u > (float)spline.getGridX1().getUmax() || // v < 0.f || v > (float)spline.getGridX2().getUmax()) { @@ -475,14 +439,14 @@ GPUdi() bool TPCFastSpaceChargeCorrection::isRealLocalInsideGrid(int32_t sector, GPUdi() void TPCFastSpaceChargeCorrection::convGridToLocal(int32_t sector, int32_t row, float gridU, float gridV, float& y, float& z) const { /// convert internal grid coordinates u,v to local y, z - getSectorRowInfo(sector, row).gridMeasured.convGridToLocal(gridU, gridV, y, z); + getRowInfo(row).gridMeasured.convGridToLocal(sector, gridU, gridV, y, z); } GPUdi() void TPCFastSpaceChargeCorrection::convRealLocalToGrid(int32_t sector, int32_t row, float y, float z, float& u, float& v, float& s) const { /// convert real y, z to the internal grid coordinates + scale - const SplineType& spline = getSpline(sector, row); - getSectorRowInfo(sector, row).gridReal.convLocalToGridUntruncated(y, z, u, v, s); + const SplineType& spline = getSplineForRow(row); + getRowInfo(row).gridReal.convLocalToGridUntruncated(sector, y, z, u, v, s); // shrink to the grid u = GPUCommonMath::Clamp(u, 0.f, (float)spline.getGridX1().getUmax()); v = GPUCommonMath::Clamp(v, 0.f, (float)spline.getGridX2().getUmax()); @@ -491,13 +455,12 @@ GPUdi() void TPCFastSpaceChargeCorrection::convRealLocalToGrid(int32_t sector, i GPUdi() void TPCFastSpaceChargeCorrection::convGridToRealLocal(int32_t sector, int32_t row, float gridU, float gridV, float& y, float& z) const { /// convert internal grid coordinates u,v to the real y, z - getSectorRowInfo(sector, row).gridReal.convGridToLocal(gridU, gridV, y, z); + getRowInfo(row).gridReal.convGridToLocal(sector, gridU, gridV, y, z); } GPUdi() void TPCFastSpaceChargeCorrection::getCorrectionLocal(int32_t sector, int32_t row, float y, float z, float& dx, float& dy, float& dz) const { - const auto& info = getSectorRowInfo(sector, row); - const SplineType& spline = getSpline(sector, row); + const SplineType& spline = getSplineForRow(row); const float* splineData = getCorrectionData(sector, row); float u, v, s; @@ -506,41 +469,38 @@ GPUdi() void TPCFastSpaceChargeCorrection::getCorrectionLocal(int32_t sector, in float dxyz[3]; spline.interpolateAtU(splineData, u, v, dxyz); - if (CAMath::Abs(dxyz[0]) > 100.f || CAMath::Abs(dxyz[1]) > 100.f || CAMath::Abs(dxyz[2]) > 100.f) { + if (CAMath::Abs(dxyz[0]) > kMaxCorrection || CAMath::Abs(dxyz[1]) > kMaxCorrection || CAMath::Abs(dxyz[2]) > kMaxCorrection) { s = 0.f; // TODO: DR: Protect from FPEs, fix upstream and remove once guaranteed that it is fixed } - dx = s * GPUCommonMath::Clamp(dxyz[0], info.minCorr[0], info.maxCorr[0]); - dy = s * GPUCommonMath::Clamp(dxyz[1], info.minCorr[1], info.maxCorr[1]); - dz = s * GPUCommonMath::Clamp(dxyz[2], info.minCorr[2], info.maxCorr[2]); + dx = s * dxyz[0]; + dy = s * dxyz[1]; + dz = s * dxyz[2]; } GPUdi() float TPCFastSpaceChargeCorrection::getCorrectionXatRealYZ(int32_t sector, int32_t row, float realY, float realZ) const { - const auto& info = getSectorRowInfo(sector, row); float u, v, s; convRealLocalToGrid(sector, row, realY, realZ, u, v, s); float dx = 0; - getSplineInvX(sector, row).interpolateAtU(getCorrectionDataInvX(sector, row), u, v, &dx); - if (CAMath::Abs(dx) > 100.f) { + getSplineInvXforRow(row).interpolateAtU(getCorrectionDataInvX(sector, row), u, v, &dx); + if (CAMath::Abs(dx) > kMaxCorrection) { s = 0.f; // TODO: DR: Protect from FPEs, fix upstream and remove once guaranteed that it is fixed } - dx = s * GPUCommonMath::Clamp(dx, info.minCorr[0], info.maxCorr[0]); - return dx; + return s * dx; } GPUdi() void TPCFastSpaceChargeCorrection::getCorrectionYZatRealYZ(int32_t sector, int32_t row, float realY, float realZ, float& y, float& z) const { float u, v, s; convRealLocalToGrid(sector, row, realY, realZ, u, v, s); - const auto& info = getSectorRowInfo(sector, row); float dyz[2]; - getSplineInvYZ(sector, row).interpolateAtU(getCorrectionDataInvYZ(sector, row), u, v, dyz); - if (CAMath::Abs(dyz[0]) > 100.f || CAMath::Abs(dyz[1]) > 100.f) { + getSplineInvYZforRow(row).interpolateAtU(getCorrectionDataInvYZ(sector, row), u, v, dyz); + if (CAMath::Abs(dyz[0]) > kMaxCorrection || CAMath::Abs(dyz[1]) > kMaxCorrection) { s = 0.f; // TODO: DR: Protect from FPEs, fix upstream and remove once guaranteed that it is fixed } - y = s * GPUCommonMath::Clamp(dyz[0], info.minCorr[1], info.maxCorr[1]); - z = s * GPUCommonMath::Clamp(dyz[1], info.minCorr[2], info.maxCorr[2]); + y = s * dyz[0]; + z = s * dyz[1]; } } // namespace o2::gpu diff --git a/GPU/TPCFastTransformation/TPCFastTransform.cxx b/GPU/TPCFastTransformation/TPCFastTransform.cxx index 2a829cfcd5471..4435505472991 100644 --- a/GPU/TPCFastTransformation/TPCFastTransform.cxx +++ b/GPU/TPCFastTransformation/TPCFastTransform.cxx @@ -113,7 +113,7 @@ void TPCFastTransform::startConstruction(const TPCFastSpaceChargeCorrection& cor mCorrection.cloneFromObject(correction, nullptr); } -void TPCFastTransform::setCalibration1(int64_t timeStamp, float t0, float vDrift) +void TPCFastTransform::setCalibration(int64_t timeStamp, float t0, float vDrift) { /// Sets all drift calibration parameters and the time stamp /// diff --git a/GPU/TPCFastTransformation/TPCFastTransform.h b/GPU/TPCFastTransformation/TPCFastTransform.h index 100c465996e7d..c8afbb57ecab8 100644 --- a/GPU/TPCFastTransformation/TPCFastTransform.h +++ b/GPU/TPCFastTransformation/TPCFastTransform.h @@ -157,7 +157,7 @@ class TPCFastTransform : public FlatObject /// /// It must be called once during construction, /// but also may be called afterwards to reset these parameters. - void setCalibration1(int64_t timeStamp, float t0, float vDrift); + void setCalibration(int64_t timeStamp, float t0, float vDrift); /// Set Lumi info void setLumi(float l) { mLumi = l; } @@ -553,7 +553,7 @@ GPUdi() void TPCFastTransform::Transform(int32_t sector, int32_t row, float pad, GPUdi() void TPCFastTransform::TransformInTimeFrame(int32_t sector, float time, float& z, float maxTimeBin) const { float l = (time - mT0 - maxTimeBin) * mVdrift; // drift length cm - z = getGeometry().convDriftLengthToZ1(sector, l); + z = getGeometry().convDriftLengthToZ(sector, l); } GPUdi() void TPCFastTransform::TransformInTimeFrame(int32_t sector, int32_t row, float pad, float time, float& x, float& y, float& z, float maxTimeBin) const @@ -592,7 +592,7 @@ GPUdi() void TPCFastTransform::TransformIdealZ(int32_t sector, float time, float /// float l = (time - mT0 - vertexTime) * mVdrift; // drift length cm - z = getGeometry().convDriftLengthToZ1(sector, l); + z = getGeometry().convDriftLengthToZ(sector, l); } GPUdi() void TPCFastTransform::TransformIdeal(int32_t sector, int32_t row, float pad, float time, float& x, float& y, float& z, float vertexTime) const diff --git a/GPU/TPCFastTransformation/TPCFastTransformGeo.h b/GPU/TPCFastTransformation/TPCFastTransformGeo.h index 23092c57b7e49..2cd145c276ea3 100644 --- a/GPU/TPCFastTransformation/TPCFastTransformGeo.h +++ b/GPU/TPCFastTransformation/TPCFastTransformGeo.h @@ -142,7 +142,6 @@ class TPCFastTransformGeo #endif GPUd() float getZmin(int32_t sector) const; GPUd() float getZmax(int32_t sector) const; - GPUd() float getZreadout(int32_t sector) const; /// _______________ Conversion of coordinate systems __________ @@ -156,10 +155,10 @@ class TPCFastTransformGeo GPUd() void convPadDriftLengthToLocal(int32_t sector, int32_t row, float pad, float driftLength, float& y, float& z) const; /// convert DriftLength -> Local c.s. - GPUd() float convDriftLengthToZ1(int32_t sector, float driftLength) const; + GPUd() float convDriftLengthToZ(int32_t sector, float driftLength) const; /// convert Z to DriftLength - GPUd() float convZtoDriftLength1(int32_t sector, float z) const; + GPUd() float convZtoDriftLength(int32_t sector, float z) const; /// convert Local c.s. -> Pad, DriftLength GPUd() void convLocalToPadDriftLength(int32_t sector, int32_t row, float y, float z, float& pad, float& l) const; @@ -262,13 +261,13 @@ GPUdi() void TPCFastTransformGeo::convPadDriftLengthToLocal(int32_t sector, int3 } } -GPUdi() float TPCFastTransformGeo::convDriftLengthToZ1(int32_t sector, float driftLength) const +GPUdi() float TPCFastTransformGeo::convDriftLengthToZ(int32_t sector, float driftLength) const { /// convert DriftLength -> Local c.s. return (sector < NumberOfSectorsA) ? (mTPCzLength - driftLength) : (driftLength - mTPCzLength); } -GPUdi() float TPCFastTransformGeo::convZtoDriftLength1(int32_t sector, float z) const +GPUdi() float TPCFastTransformGeo::convZtoDriftLength(int32_t sector, float z) const { /// convert Z to DriftLength return (sector < NumberOfSectorsA) ? (mTPCzLength - z) : (z + mTPCzLength); @@ -294,16 +293,6 @@ GPUdi() float TPCFastTransformGeo::getZmax(int32_t sector) const } } -GPUdi() float TPCFastTransformGeo::getZreadout(int32_t sector) const -{ - /// z readout for the sector - if (sector < NumberOfSectorsA) { // TPC side A - return mTPCzLength; - } else { // TPC side C - return -mTPCzLength; - } -} - GPUdi() void TPCFastTransformGeo::convLocalToPadDriftLength(int32_t sector, int32_t row, float y, float z, float& pad, float& l) const { /// convert Local c.s. -> Pad, DriftLength diff --git a/GPU/TPCFastTransformation/TPCFastTransformGeoPOD.h b/GPU/TPCFastTransformation/TPCFastTransformGeoPOD.h index dbb6176dd47b9..3cbf1dae0c4dd 100644 --- a/GPU/TPCFastTransformation/TPCFastTransformGeoPOD.h +++ b/GPU/TPCFastTransformation/TPCFastTransformGeoPOD.h @@ -91,13 +91,13 @@ struct TPCFastTransformGeoPOD { } /// convert DriftLength -> Local c.s. - inline static constexpr float convDriftLengthToZ1(uint32_t sector, float driftLength) + inline static constexpr float convDriftLengthToZ(uint32_t sector, float driftLength) { return (sector < getNumberOfSectorsA()) ? (getTPCzLength() - driftLength) : (driftLength - getTPCzLength()); } /// convert Z to DriftLength - inline static constexpr float convZtoDriftLength1(uint32_t sector, float z) + inline static constexpr float convZtoDriftLength(uint32_t sector, float z) { return (sector < getNumberOfSectorsA()) ? (getTPCzLength() - z) : (z + getTPCzLength()); } diff --git a/GPU/TPCFastTransformation/TPCFastTransformPOD.cxx b/GPU/TPCFastTransformation/TPCFastTransformPOD.cxx index 58635995a99c1..9c8f715da94ec 100644 --- a/GPU/TPCFastTransformation/TPCFastTransformPOD.cxx +++ b/GPU/TPCFastTransformation/TPCFastTransformPOD.cxx @@ -69,7 +69,7 @@ size_t TPCFastTransformPOD::estimateSize(const TPCFastSpaceChargeCorrection& ori for (int is = 0; is < 3; is++) { for (int sector = 0; sector < origCorr.mGeo.getNumberOfSectors(); sector++) { for (int row = 0; row < NROWS; row++) { - const auto& spline = origCorr.getSpline(sector, row); + const auto& spline = origCorr.getSplineForRow(row); int nPar = spline.getNumberOfParameters(); if (is == 1) { nPar = nPar / 3; @@ -96,7 +96,7 @@ void TPCFastTransformPOD::print() const } } const size_t scenOffset = getScenarioOffset(0); - const auto& spline = getSpline(0, 0); + const auto& spline = getSplineForRow(0); LOGP(info, "scenOffset={} spline_addr={:p} expected={:p}", scenOffset, (void*)&spline, (void*)(getThis() + scenOffset)); const float* splineData = getCorrectionData(0, 0); @@ -120,11 +120,11 @@ TPCFastTransformPOD* TPCFastTransformPOD::create(char* buff, size_t buffSize, co // copy fixed size data --- start podMap.mNumberOfScenarios = origCorr.mNumberOfScenarios; - for (int sector = 0; sector < TPCFastTransformGeo::getNumberOfSectors(); sector++) { - for (int row = 0; row < NROWS; row++) { - podMap.mSectorRowInfos[NROWS * sector + row] = origCorr.getSectorRowInfo(sector, row); - } + + for (int row = 0; row < NROWS; row++) { + podMap.mRowInfos[row] = origCorr.getRowInfo(row); } + podMap.mTimeStamp = origCorr.mTimeStamp; // // init data members coming from the TPCFastTrasform @@ -169,37 +169,27 @@ TPCFastTransformPOD* TPCFastTransformPOD::create(char* buff, size_t buffSize, co nextDynOffs = alignOffset(nextDynOffs + spline.getFlatBufferSize()); } - // copy splines data + // copy spline data for (int is = 0; is < 3; is++) { float* data = reinterpret_cast(buff + nextDynOffs); LOGP(debug, "splinID={} start offset {} -> {}", is, nextDynOffs, (void*)data); + + // metadata + size_t sectorDataSizeBytes = origCorr.mSectorDataSizeBytes[is]; + for (int sector = 0; sector < origCorr.mGeo.getNumberOfSectors(); sector++) { - podMap.mSplineDataOffsets[sector][is] = nextDynOffs; - size_t rowDataOffs = 0; - for (int row = 0; row < NROWS; row++) { - const auto& spline = origCorr.getSpline(sector, row); - const float* dataOr = origCorr.getCorrectionData(sector, 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 sector:{} row:{} to offset {}", nPar, is, sector, 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 sector{}/row{} to {}, overflowing the buffer of size {}", nbcopy, is, sector, row, nextDynOffs, buffSize)); - } - std::memcpy(data, dataOr, nbcopy); - podMap.getSectorRowInfo(sector, row).dataOffsetBytes[is] = rowDataOffs; - rowDataOffs += nbcopy; - data += nPar; - nextDynOffs += nbcopy; - } + podMap.mSplineDataOffsets[sector][is] = sectorDataSizeBytes * sector; } + if (buffSize < nextDynOffs + sectorDataSizeBytes * origCorr.mGeo.getNumberOfSectors()) { + throw std::runtime_error(fmt::format("attempt to copy {} bytes of data for spline{} to {}, overflowing the buffer of size {}", sectorDataSizeBytes, is, nextDynOffs, buffSize)); + } + const char* dataOr = origCorr.mCorrectionData[is]; + size_t dataSize = origCorr.mGeo.getNumberOfSectors() * sectorDataSizeBytes; + std::memcpy(data, dataOr, dataSize); + nextDynOffs += dataSize; } - podMap.mTotalSize = alignOffset(nextDynOffs); + + podMap.mTotalSize = nextDynOffs; if (buffSize != podMap.mTotalSize) { throw std::runtime_error(fmt::format("Estimated buffer size {} differs from filled one {}", buffSize, podMap.mTotalSize)); } diff --git a/GPU/TPCFastTransformation/TPCFastTransformPOD.h b/GPU/TPCFastTransformation/TPCFastTransformPOD.h index b30a3d52ec696..c7e06d4b47ca4 100644 --- a/GPU/TPCFastTransformation/TPCFastTransformPOD.h +++ b/GPU/TPCFastTransformation/TPCFastTransformPOD.h @@ -52,7 +52,7 @@ class TPCFastTransformPOD public: using SliceInfo = TPCFastSpaceChargeCorrection::SliceInfo; // obsolete using GridInfo = TPCFastSpaceChargeCorrection::GridInfo; - using SectorRowInfo = TPCFastSpaceChargeCorrection::SectorRowInfo; + using RowInfo = TPCFastSpaceChargeCorrection::RowInfo; using SplineTypeXYZ = TPCFastSpaceChargeCorrection::SlimSplineTypeXYZ; using SplineTypeInvX = TPCFastSpaceChargeCorrection::SlimSplineTypeInvX; @@ -112,10 +112,10 @@ class TPCFastTransformPOD GPUd() const TPCFastTransformGeoPOD& getGeometry() const { return mGeo; } /// Gives TPC sector & row info - GPUd() const SectorRowInfo& getSectorRowInfo(int32_t sector, int32_t row) const { return mSectorRowInfos[NROWS * sector + row]; } + GPUd() const RowInfo& getRowInfo(int32_t row) const { return mRowInfos[row]; } /// Gives TPC sector & row info - GPUd() SectorRowInfo& getSectorRowInfo(int32_t sector, int32_t row) { return mSectorRowInfos[NROWS * sector + row]; } + GPUd() RowInfo& getRowInfo(int32_t row) { return mRowInfos[row]; } /// Gives its own size including dynamic part GPUd() size_t size() const { return mTotalSize; } @@ -159,22 +159,22 @@ class TPCFastTransformPOD /// Sets CTP Lumi estimator GPUd() void setLumi(float v) { mLumi = v; } - GPUd() void setCalibration1(int64_t timeStamp, float t0, float vDrift); + GPUd() void setCalibration(int64_t timeStamp, float t0, float vDrift); /// Gives a reference to a spline - GPUd() const SplineType& getSpline(int32_t sector, int32_t row) const { return *reinterpret_cast(getThis() + getScenarioOffset(getSectorRowInfo(sector, row).splineScenarioID)); } + GPUd() const SplineType& getSplineForRow(int32_t row) const { return *reinterpret_cast(getThis() + getScenarioOffset(getRowInfo(row).splineScenarioID)); } /// Gives pointer to spline data - GPUd() const float* getCorrectionData(int32_t sector, int32_t row, int32_t iSpline = 0) const { return reinterpret_cast(getThis() + mSplineDataOffsets[sector][iSpline] + getSectorRowInfo(sector, row).dataOffsetBytes[iSpline]); } + GPUd() const float* getCorrectionData(int32_t sector, int32_t row, int32_t iSpline = 0) const { return reinterpret_cast(getThis() + mSplineDataOffsets[sector][iSpline] + getRowInfo(row).dataOffsetBytes[iSpline]); } /// Gives const pointer to a spline for the inverse X correction - GPUd() const SplineTypeInvX& getSplineInvX(int32_t sector, int32_t row) const { return reinterpret_cast(getSpline(sector, row)); } + GPUd() const SplineTypeInvX& getSplineInvXforRow(int32_t row) const { return reinterpret_cast(getSplineForRow(row)); } /// Gives pointer to spline data for the inverse X correction GPUd() const float* getCorrectionDataInvX(int32_t sector, int32_t row) const { return getCorrectionData(sector, row, 1); } /// Gives const pointer to a spline for the inverse YZ correction - GPUd() const SplineTypeInvYZ& getSplineInvYZ(int32_t sector, int32_t row) const { return reinterpret_cast(getSpline(sector, row)); } + GPUd() const SplineTypeInvYZ& getSplineInvYZforRow(int32_t row) const { return reinterpret_cast(getSplineForRow(row)); } /// Gives pointer to spline data for the inverse YZ correction GPUd() const float* getCorrectionDataInvYZ(int32_t sector, int32_t row) const { return getCorrectionData(sector, row, 2); } @@ -242,6 +242,7 @@ class TPCFastTransformPOD static constexpr int NROWS = o2::tpc::constants::MAXGLOBALPADROW; static constexpr int NSECTORS = o2::tpc::constants::MAXSECTOR; + static constexpr int NSECTORSA = o2::tpc::constants::MAXSECTOR / 2; static constexpr int NSplineIDs = 3; ///< number of spline data sets for each sector/row private: @@ -295,16 +296,16 @@ class TPCFastTransformPOD float mIDC; ///< IDC estimator (for info only) TPCFastTransformGeoPOD mGeo; ///< TPC geometry information - SectorRowInfo mSectorRowInfos[NROWS * TPCFastTransformGeo::getNumberOfSectors()]; + RowInfo mRowInfos[NROWS]; ClassDefNV(TPCFastTransformPOD, 0); }; GPUdi() void TPCFastTransformPOD::getCorrectionLocal(int32_t sector, int32_t row, float y, float z, float& dx, float& dy, float& dz) const { - const auto& info = getSectorRowInfo(sector, row); + const auto& info = getRowInfo(row); const int32_t isc = info.splineScenarioID; - const SplineType& spline = getSpline(sector, row); + const SplineType& spline = getSplineForRow(row); const float* splineData = getCorrectionData(sector, row); float u, v, s; @@ -316,32 +317,32 @@ GPUdi() void TPCFastTransformPOD::getCorrectionLocal(int32_t sector, int32_t row float dxyz[3]; spline.interpolateAtUZeroCopy(g1buf, g2buf, splineData, u, v, dxyz); - if (CAMath::Abs(dxyz[0]) > 100.f || CAMath::Abs(dxyz[1]) > 100.f || CAMath::Abs(dxyz[2]) > 100.f) { + if (CAMath::Abs(dxyz[0]) > TPCFastSpaceChargeCorrection::kMaxCorrection || CAMath::Abs(dxyz[1]) > TPCFastSpaceChargeCorrection::kMaxCorrection || CAMath::Abs(dxyz[2]) > TPCFastSpaceChargeCorrection::kMaxCorrection) { s = 0.f; // TODO: DR: Protect from FPEs, fix upstream and remove once guaranteed that it is fixed } - dx = s * GPUCommonMath::Clamp(dxyz[0], info.minCorr[0], info.maxCorr[0]); - dy = s * GPUCommonMath::Clamp(dxyz[1], info.minCorr[1], info.maxCorr[1]); - dz = s * GPUCommonMath::Clamp(dxyz[2], info.minCorr[2], info.maxCorr[2]); + dx = s * dxyz[0]; + dy = s * dxyz[1]; + dz = s * dxyz[2]; } GPUdi() float TPCFastTransformPOD::getCorrectionXatRealYZ(int32_t sector, int32_t row, float realY, float realZ) const { - const auto& info = getSectorRowInfo(sector, row); + const auto& info = getRowInfo(row); float u, v, s; convRealLocalToGrid(sector, row, realY, realZ, u, v, s); const int32_t isc = info.splineScenarioID; - const auto& spline = getSplineInvX(sector, row); + const auto& spline = getSplineInvXforRow(row); const char* g1buf = getSplineFlatBuffer(isc); const char* g2buf = getGridX2FlatBuffer(spline, isc); float dx = 0; spline.interpolateAtUZeroCopy(g1buf, g2buf, getCorrectionDataInvX(sector, row), u, v, &dx); - if (CAMath::Abs(dx) > 100.f) { + if (CAMath::Abs(dx) > TPCFastSpaceChargeCorrection::kMaxCorrection) { s = 0.f; // TODO: DR: Protect from FPEs, fix upstream and remove once guaranteed that it is fixed } - dx = s * GPUCommonMath::Clamp(dx, info.minCorr[0], info.maxCorr[0]); + dx = s * dx; return dx; } @@ -349,27 +350,27 @@ GPUdi() void TPCFastTransformPOD::getCorrectionYZatRealYZ(int32_t sector, int32_ { float u, v, s; convRealLocalToGrid(sector, row, realY, realZ, u, v, s); - const auto& info = getSectorRowInfo(sector, row); + const auto& info = getRowInfo(row); const int32_t isc = info.splineScenarioID; - const auto& spline = getSplineInvYZ(sector, row); + const auto& spline = getSplineInvYZforRow(row); const char* g1buf = getSplineFlatBuffer(isc); const char* g2buf = getGridX2FlatBuffer(spline, isc); float dyz[2]; spline.interpolateAtUZeroCopy(g1buf, g2buf, getCorrectionDataInvYZ(sector, row), u, v, dyz); - if (CAMath::Abs(dyz[0]) > 100.f || CAMath::Abs(dyz[1]) > 100.f) { + if (CAMath::Abs(dyz[0]) > TPCFastSpaceChargeCorrection::kMaxCorrection || CAMath::Abs(dyz[1]) > TPCFastSpaceChargeCorrection::kMaxCorrection) { s = 0.f; // TODO: DR: Protect from FPEs, fix upstream and remove once guaranteed that it is fixed } - y = s * GPUCommonMath::Clamp(dyz[0], info.minCorr[1], info.maxCorr[1]); - z = s * GPUCommonMath::Clamp(dyz[1], info.minCorr[2], info.maxCorr[2]); + y = s * dyz[0]; + z = s * dyz[1]; } GPUdi() void TPCFastTransformPOD::convLocalToGrid(int32_t sector, int32_t row, float y, float z, float& u, float& v, float& s) const { /// convert local y, z to internal grid coordinates u,v /// return values: u, v, scaling factor - const SplineType& spline = getSpline(sector, row); - getSectorRowInfo(sector, row).gridMeasured.convLocalToGridUntruncated(y, z, u, v, s); + const SplineType& spline = getSplineForRow(row); + getRowInfo(row).gridMeasured.convLocalToGridUntruncated(sector, y, z, u, v, s); // shrink to the grid u = GPUCommonMath::Clamp(u, 0.f, (float)spline.getGridX1().getUmax()); v = GPUCommonMath::Clamp(v, 0.f, (float)spline.getGridX2().getUmax()); @@ -378,14 +379,14 @@ GPUdi() void TPCFastTransformPOD::convLocalToGrid(int32_t sector, int32_t row, f GPUdi() void TPCFastTransformPOD::convGridToLocal(int32_t sector, int32_t row, float gridU, float gridV, float& y, float& z) const { /// convert internal grid coordinates u,v to local y, z - getSectorRowInfo(sector, row).gridMeasured.convGridToLocal(gridU, gridV, y, z); + getRowInfo(row).gridMeasured.convGridToLocal(sector, gridU, gridV, y, z); } GPUdi() void TPCFastTransformPOD::convRealLocalToGrid(int32_t sector, int32_t row, float y, float z, float& u, float& v, float& s) const { /// convert real y, z to the internal grid coordinates + scale - const SplineType& spline = getSpline(sector, row); - getSectorRowInfo(sector, row).gridReal.convLocalToGridUntruncated(y, z, u, v, s); + const SplineType& spline = getSplineForRow(row); + getRowInfo(row).gridReal.convLocalToGridUntruncated(sector, y, z, u, v, s); // shrink to the grid u = GPUCommonMath::Clamp(u, 0.f, (float)spline.getGridX1().getUmax()); v = GPUCommonMath::Clamp(v, 0.f, (float)spline.getGridX2().getUmax()); @@ -394,15 +395,15 @@ GPUdi() void TPCFastTransformPOD::convRealLocalToGrid(int32_t sector, int32_t ro GPUdi() void TPCFastTransformPOD::convGridToRealLocal(int32_t sector, int32_t row, float gridU, float gridV, float& y, float& z) const { /// convert internal grid coordinates u,v to the real y, z - getSectorRowInfo(sector, row).gridReal.convGridToLocal(gridU, gridV, y, z); + getRowInfo(row).gridReal.convGridToLocal(sector, gridU, gridV, y, z); } GPUdi() bool TPCFastTransformPOD::isLocalInsideGrid(int32_t sector, int32_t row, float y, float z) const { /// check if local y, z are inside the grid float u, v, s; - getSectorRowInfo(sector, row).gridMeasured.convLocalToGridUntruncated(y, z, u, v, s); - const auto& spline = getSpline(sector, row); + getRowInfo(row).gridMeasured.convLocalToGridUntruncated(sector, y, z, u, v, s); + const auto& spline = getSplineForRow(row); // shrink to the grid if (u < 0.f || u > (float)spline.getGridX1().getUmax() || // v < 0.f || v > (float)spline.getGridX2().getUmax()) { @@ -415,8 +416,8 @@ GPUdi() bool TPCFastTransformPOD::isRealLocalInsideGrid(int32_t sector, int32_t { /// check if local y, z are inside the grid float u, v, s; - getSectorRowInfo(sector, row).gridReal.convLocalToGridUntruncated(y, z, u, v, s); - const auto& spline = getSpline(sector, row); + getRowInfo(row).gridReal.convLocalToGridUntruncated(sector, y, z, u, v, s); + const auto& spline = getSplineForRow(row); // shrink to the grid if (u < 0.f || u > (float)spline.getGridX1().getUmax() || // v < 0.f || v > (float)spline.getGridX2().getUmax()) { @@ -500,7 +501,7 @@ GPUdi() void TPCFastTransformPOD::TransformXYZ(int32_t sector, int32_t row, floa GPUdi() void TPCFastTransformPOD::TransformInTimeFrame(int32_t sector, float time, float& z, float maxTimeBin) const { float l = (time - mT0 - maxTimeBin) * mVdrift; // drift length cm - z = getGeometry().convDriftLengthToZ1(sector, l); + z = getGeometry().convDriftLengthToZ(sector, l); } GPUdi() void TPCFastTransformPOD::TransformInTimeFrame(int32_t sector, int32_t row, float pad, float time, float& x, float& y, float& z, float maxTimeBin) const @@ -538,7 +539,7 @@ GPUdi() void TPCFastTransformPOD::TransformIdealZ(int32_t sector, float time, fl /// float l = (time - mT0 - vertexTime) * mVdrift; // drift length cm - z = getGeometry().convDriftLengthToZ1(sector, l); + z = getGeometry().convDriftLengthToZ(sector, l); } GPUdi() void TPCFastTransformPOD::TransformIdeal(int32_t sector, int32_t row, float pad, float time, float& x, float& y, float& z, float vertexTime) const @@ -692,7 +693,7 @@ GPUdi() float TPCFastTransformPOD::convVertexTimeToZOffset(int32_t sector, float } #ifndef GPUCA_GPUCODE_DEVICE // Functions not needed during GPU processing -GPUdi() void TPCFastTransformPOD::setCalibration1(int64_t timeStamp, float t0, float vDrift) +GPUdi() void TPCFastTransformPOD::setCalibration(int64_t timeStamp, float t0, float vDrift) { mTimeStamp = timeStamp; mT0 = t0; diff --git a/GPU/TPCFastTransformation/TPCFastTransformationLinkDef_O2.h b/GPU/TPCFastTransformation/TPCFastTransformationLinkDef_O2.h index cc12badb1e654..c9f89f3dcac76 100644 --- a/GPU/TPCFastTransformation/TPCFastTransformationLinkDef_O2.h +++ b/GPU/TPCFastTransformation/TPCFastTransformationLinkDef_O2.h @@ -47,10 +47,10 @@ #pragma link C++ class o2::gpu::TPCFastSpaceChargeCorrection + ; #pragma link C++ class o2::gpu::TPCFastSpaceChargeCorrection::SliceInfo + ; -#pragma link C++ class o2::gpu::TPCFastSpaceChargeCorrection::SectorRowInfo + ; +#pragma link C++ class o2::gpu::TPCFastSpaceChargeCorrection::RowInfo + ; #pragma link C++ class o2::gpu::TPCFastSpaceChargeCorrection::GridInfo + ; #pragma read sourceClass = "o2::gpu::TPCFastSpaceChargeCorrection" targetClass = "o2::gpu::TPCFastSpaceChargeCorrection" source = "o2::gpu::TPCFastSpaceChargeCorrection::SliceInfo mSliceInfo[36]" version = "[-3]" target = "" code = "{}"; -#pragma read sourceClass = "o2::gpu::TPCFastSpaceChargeCorrection" targetClass = "o2::gpu::TPCFastSpaceChargeCorrection" source = "size_t mSliceDataSizeBytes[3]" version = "[-3]" target = "mCorrectionDataSize" code = "{ for (int i=0; i<3; i++) mCorrectionDataSize[i] = onfile.mSliceDataSizeBytes[i] * o2::tpc::constants::MAXSECTOR; }"; +#pragma read sourceClass = "o2::gpu::TPCFastSpaceChargeCorrection" targetClass = "o2::gpu::TPCFastSpaceChargeCorrection" source = "size_t mSliceDataSizeBytes[3]" version = "[-3]" target = "mSectorDataSizeBytes" code = "{ for (int i=0; i<3; i++) mSectorDataSizeBytes[i] = onfile.mSliceDataSizeBytes[i]; }"; #pragma read sourceClass = "o2::gpu::TPCFastSpaceChargeCorrection" targetClass = "o2::gpu::TPCFastSpaceChargeCorrection" source = "float fInterpolationSafetyMargin" version = "[-3]" target = "" code = "{}"; #pragma link C++ struct o2::gpu::TPCSlowSpaceChargeCorrection + ; diff --git a/GPU/TPCFastTransformation/macro/TPCFastTransformInit.C b/GPU/TPCFastTransformation/macro/TPCFastTransformInit.C index bc6fafbaa8bd0..dd8bb1a789db9 100644 --- a/GPU/TPCFastTransformation/macro/TPCFastTransformInit.C +++ b/GPU/TPCFastTransformation/macro/TPCFastTransformInit.C @@ -180,7 +180,8 @@ void TPCFastTransformInit(const char* fileName = "debugVoxRes.root", const char* const char* fileName = outFileName; - // fileName = "~/test/master/TPCFastTransform_VoxRes.root"; + // file with the old data format + // fileName = "~/alidock/test/master/TPCFastTransform_VoxRes.root"; std::cout << "load corrections from file " << fileName << std::endl; @@ -193,17 +194,14 @@ void TPCFastTransformInit(const char* fileName = "debugVoxRes.root", const char* const o2::gpu::TPCFastTransformGeo& geo = helper->getGeometry(); - // for (int32_t iSector = 0; iSector < geo.getNumberOfSectors(); iSector++) { - for (int32_t iSector = 0; iSector < 1; iSector++) { - for (int32_t iRow = 0; iRow < geo.getNumberOfRows(); iRow++) { - auto& info = corr.getSectorRowInfo(iSector, iRow); - std::cout << "sector " << iSector << " row " << iRow - << " gridY0 " << info.gridMeasured.getY0() << " gridZ0 " << info.gridMeasured.getZ0() - << " scaleYtoGrid " << info.gridMeasured.getYscale() << " scaleLtoGrid " << info.gridMeasured.getZscale() - << " gridRealY0 " << info.gridReal.getY0() << " gridRealZ0 " << info.gridReal.getZ0() - << " scaleRealYtoGrid " << info.gridReal.getYscale() << " scaleRealLtoGrid " << info.gridReal.getZscale() - << std::endl; - } + for (int32_t iRow = 0; iRow < geo.getNumberOfRows(); iRow++) { + auto& info = corr.getRowInfo(iRow); + std::cout << " row " << iRow + << " gridY0 " << info.gridMeasured.getY0() << " gridZ0 " << info.gridMeasured.getZ0() + << " scaleYtoGrid " << info.gridMeasured.getYscale() << " scaleLtoGrid " << info.gridMeasured.getZscale() + << " gridRealY0 " << info.gridReal.getY0() << " gridRealZ0 " << info.gridReal.getZ0() + << " scaleRealYtoGrid " << info.gridReal.getYscale() << " scaleRealLtoGrid " << info.gridReal.getZscale() + << std::endl; } } } @@ -307,17 +305,12 @@ void TPCFastTransformInit(const char* fileName = "debugVoxRes.root", const char* auto getInvCorrections = [&](int iSector, int iRow, float realY, float realZ, float& ix, float& iy, float& iz) { // get the inverse corrections ix, iy, iz at x,y,z ix = corr.getCorrectionXatRealYZ(iSector, iRow, realY, realZ); - const auto c = corr.getCorrectionYZatRealYZ(iSector, iRow, realY, realZ); - iy = c[0]; - iz = c[1]; + corr.getCorrectionYZatRealYZ(iSector, iRow, realY, realZ, iy, iz); }; auto getAllCorrections = [&](int iSector, int iRow, float y, float z, float& cx, float& cy, float& cz, float& ix, float& iy, float& iz) { // get the corrections cx,cy,cz at x,y,z - const auto c = corr.getCorrectionLocal(iSector, iRow, y, z); - cx = c[0]; - cy = c[1]; - cz = c[2]; + corr.getCorrectionLocal(iSector, iRow, y, z, cx, cy, cz); getInvCorrections(iSector, iRow, y + cy, z + cz, ix, iy, iz); }; @@ -477,10 +470,10 @@ void TPCFastTransformInit(const char* fileName = "debugVoxRes.root", const char* // the spline grid - const auto& gridY = corr.getSpline(iSector, iRow).getGridX1(); - const auto& gridZ = corr.getSpline(iSector, iRow).getGridX2(); + const auto& gridY = corr.getSplineForRow(iRow).getGridX1(); + const auto& gridZ = corr.getSplineForRow(iRow).getGridX2(); if (iSector == 0 && iRow == 0) { - std::cout << "spline scenario " << corr.getSectorRowInfo(iSector, iRow).splineScenarioID << std::endl; + std::cout << "spline scenario " << corr.getRowInfo(iRow).splineScenarioID << std::endl; std::cout << "spline grid Y: u = " << 0 << ".." << gridY.getUmax() << ", x = " << gridY.getXmin() << ".." << gridY.getXmax() << std::endl; std::cout << "spline grid Z: u = " << 0 << ".." << gridZ.getUmax() << ", x = " << gridZ.getXmin() << ".." << gridZ.getXmax() << std::endl; } @@ -582,15 +575,15 @@ void TPCFastTransformInit(const char* fileName = "debugVoxRes.root", const char* float correctionY = point.mDy; float correctionZ = point.mDz; if (direction == 0) { - auto [cx, cy, cz] = - corr.getCorrectionLocal(iSector, iRow, y, z); + float cx, cy, cz; + corr.getCorrectionLocal(iSector, iRow, y, z, cx, cy, cz); ntFitPoints->Fill(iSector, iRow, x, y, z, correctionX, correctionY, correctionZ, cx, cy, cz); } else { float cx = corr.getCorrectionXatRealYZ(iSector, iRow, y, z); - auto [cy, cz] = - corr.getCorrectionYZatRealYZ(iSector, iRow, y, z); + float cy, cz; + corr.getCorrectionYZatRealYZ(iSector, iRow, y, z, cy, cz); ntInvFitPoints->Fill(iSector, iRow, x, y, z, correctionX, correctionY, correctionZ, cx, cy, cz); }