Skip to content

Commit 33f4a4c

Browse files
ekryshenekryshen
committed
Fixes for modules at negative eta
1 parent 879a535 commit 33f4a4c

File tree

4 files changed

+86
-38
lines changed

4 files changed

+86
-38
lines changed

‎Detectors/Upgrades/ALICE3/ECal/base/include/ECalBase/Geometry.h‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -61,6 +61,7 @@ class Geometry
6161
double getSamplingAlpha() { return mSamplingAlpha; }
6262
double getCrystalDeltaPhi() { return 2 * std::atan(mCrystalModW / 2 / mRMin); }
6363
double getSamplingDeltaPhi() { return 2 * std::atan(mSamplingModW / 2 / mRMin); }
64+
double getFrontFaceMaxEta(int i);
6465
double getCrystalPhiMin();
6566
double getSamplingPhiMin();
6667
int getNModulesZ() { return mNModulesZ; }

‎Detectors/Upgrades/ALICE3/ECal/base/src/Geometry.cxx‎

Lines changed: 11 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -73,6 +73,12 @@ double Geometry::getSamplingPhiMin()
7373
return (superModuleDeltaPhi - samplingDeltaPhi * mNSamplingModulesPhi) / 2.;
7474
}
7575

76+
double Geometry::getFrontFaceMaxEta(int i)
77+
{
78+
double theta = std::atan(mRMin / getFrontFaceZatMinR(i));
79+
return -std::log(std::tan(theta / 2.));
80+
}
81+
7682
//==============================================================================
7783
void Geometry::fillFrontFaceCenterCoordinates()
7884
{
@@ -153,7 +159,7 @@ int Geometry::getCellID(int moduleId, int sectorId, bool isCrystal)
153159
if (sectorId % 2 == 0) { // sampling at positive eta
154160
cellID = sectorId / 2 * mNModulesZ + moduleId + mNSamplingModulesZ + mNCrystalModulesZ * 2;
155161
} else { // sampling at negative eta
156-
cellID = sectorId / 2 * mNModulesZ - moduleId + mNSamplingModulesZ;
162+
cellID = sectorId / 2 * mNModulesZ - moduleId + mNSamplingModulesZ - 1;
157163
}
158164
}
159165
return cellID;
@@ -206,13 +212,15 @@ void Geometry::detIdToGlobalPosition(int detId, double& x, double& y, double& z)
206212
{
207213
int chamber, sector, iphi, iz;
208214
detIdToRelIndex(detId, chamber, sector, iphi, iz);
215+
double r = 0;
209216
if (iz < mNSamplingModulesZ + mNCrystalModulesZ) {
210217
z = -mFrontFaceCenterZ[mNSamplingModulesZ + mNCrystalModulesZ - iz - 1];
218+
r = mFrontFaceCenterR[mNSamplingModulesZ + mNCrystalModulesZ - iz - 1];
211219
} else {
212-
z = +mFrontFaceCenterZ[iz % (mNSamplingModulesZ + mNCrystalModulesZ)];
220+
z = mFrontFaceCenterZ[iz % (mNSamplingModulesZ + mNCrystalModulesZ)];
221+
r = mFrontFaceCenterR[iz % (mNSamplingModulesZ + mNCrystalModulesZ)];
213222
}
214223
double phi = chamber == 1 ? mFrontFaceCenterCrystalPhi[iphi] : mFrontFaceCenterSamplingPhi[iphi];
215-
double r = mFrontFaceCenterR[iz % (mNSamplingModulesZ + mNCrystalModulesZ)];
216224
x = r * std::cos(phi);
217225
y = r * std::sin(phi);
218226
}

‎Detectors/Upgrades/ALICE3/ECal/reconstruction/include/ECalReconstruction/Clusterizer.h‎

Lines changed: 25 additions & 17 deletions
Original file line numberDiff line numberDiff line change
@@ -48,25 +48,33 @@ class Clusterizer
4848
void setClusteringThreshold(double threshold) { mClusteringThreshold = threshold; }
4949
void setCrystalDigitThreshold(double threshold) { mCrystalDigitThreshold = threshold; }
5050
void setSamplingDigitThreshold(double threshold) { mSamplingDigitThreshold = threshold; }
51+
void setCrystalEnergyCorrectionPars(std::vector<double> pars) { mCrystalEnergyCorrectionPars = pars; }
52+
void setSamplingEnergyCorrectionPars(std::vector<double> pars) { mSamplingEnergyCorrectionPars = pars; }
53+
void setCrystalZCorrectionPars(std::vector<double> pars) { mCrystalZCorrectionPars = pars; }
54+
void setSamplingZCorrectionPars(std::vector<double> pars) { mSamplingZCorrectionPars = pars; }
5155

5256
private:
53-
std::vector<std::vector<int>> mDigitIndices; // 2D map of digit indices used for recursive cluster finding
54-
bool mUnfoldClusters = true; // to perform cluster unfolding
55-
double mCrystalDigitThreshold = 0.040; // minimal energy of crystal digit
56-
double mSamplingDigitThreshold = 0.100; // minimal energy of sampling digit
57-
double mClusteringThreshold = 0.050; // minimal energy of digit to start clustering (GeV)
58-
double mClusteringTimeGate = 1e9; // maximal time difference between digits to be accepted to clusters (in ns)
59-
int mNLMMax = 30; // maximal number of local maxima in unfolding
60-
double mLogWeight = 4.; // cutoff used in log. weight calculation
61-
double mUnfogingEAccuracy = 1.e-4; // accuracy of energy calculation in unfoding prosedure (GeV)
62-
double mUnfogingXZAccuracy = 1.e-2; // accuracy of position calculation in unfolding procedure (cm)
63-
int mNMaxIterations = 100; // maximal number of iterations in unfolding procedure
64-
double mLocalMaximumCut = 0.015; // minimal height of local maximum over neighbours
65-
bool mApplyCorrectionZ = 1; // z-correction
66-
bool mApplyCorrectionE = 1; // energy-correction
67-
TF1* fCrystalShowerShape; //! Crystal shower shape
68-
TF1* fSamplingShowerShape; //! Sampling shower shape
69-
TF1* fCrystalRMS; //! Crystal RMS
57+
std::vector<std::vector<int>> mDigitIndices; // 2D map of digit indices used for recursive cluster finding
58+
bool mUnfoldClusters = true; // to perform cluster unfolding
59+
double mCrystalDigitThreshold = 0.040; // minimal energy of crystal digit
60+
double mSamplingDigitThreshold = 0.100; // minimal energy of sampling digit
61+
double mClusteringThreshold = 0.050; // minimal energy of digit to start clustering (GeV)
62+
double mClusteringTimeGate = 1e9; // maximal time difference between digits to be accepted to clusters (in ns)
63+
int mNLMMax = 30; // maximal number of local maxima in unfolding
64+
double mLogWeight = 4.; // cutoff used in log. weight calculation
65+
double mUnfogingEAccuracy = 1.e-4; // accuracy of energy calculation in unfoding prosedure (GeV)
66+
double mUnfogingXZAccuracy = 1.e-2; // accuracy of position calculation in unfolding procedure (cm)
67+
int mNMaxIterations = 100; // maximal number of iterations in unfolding procedure
68+
double mLocalMaximumCut = 0.015; // minimal height of local maximum over neighbours
69+
bool mApplyCorrectionZ = 1; // apply z-correction
70+
bool mApplyCorrectionE = 1; // apply energy-correction
71+
TF1* fCrystalShowerShape; //! Crystal shower shape
72+
TF1* fSamplingShowerShape; //! Sampling shower shape
73+
TF1* fCrystalRMS; //! Crystal RMS
74+
std::vector<double> mCrystalEnergyCorrectionPars; // crystal energy-correction parameters
75+
std::vector<double> mSamplingEnergyCorrectionPars; // sampling energy-correction parameters
76+
std::vector<double> mCrystalZCorrectionPars; // crystal z-correction parameters
77+
std::vector<double> mSamplingZCorrectionPars; // sampling z-correction parameters
7078
};
7179

7280
} // namespace ecal

‎Detectors/Upgrades/ALICE3/ECal/reconstruction/src/Clusterizer.cxx‎

Lines changed: 49 additions & 18 deletions
Original file line numberDiff line numberDiff line change
@@ -31,6 +31,45 @@ Clusterizer::Clusterizer(bool applyCorrectionZ, bool applyCorrectionE)
3131
mDigitIndices.resize(geo.getNrows(), std::vector<int>(geo.getNcols(), -1));
3232
mApplyCorrectionZ = applyCorrectionZ;
3333
mApplyCorrectionE = applyCorrectionE;
34+
35+
mCrystalEnergyCorrectionPars.reserve(6);
36+
mCrystalEnergyCorrectionPars[0] = 0.00444;
37+
mCrystalEnergyCorrectionPars[1] = -1.322;
38+
mCrystalEnergyCorrectionPars[2] = 1.021;
39+
mCrystalEnergyCorrectionPars[3] = 0.0018;
40+
mCrystalEnergyCorrectionPars[4] = 0.;
41+
mCrystalEnergyCorrectionPars[5] = 0.;
42+
43+
mSamplingEnergyCorrectionPars.reserve(6);
44+
mSamplingEnergyCorrectionPars[0] = 0.0033;
45+
mSamplingEnergyCorrectionPars[1] = -2.09;
46+
mSamplingEnergyCorrectionPars[2] = 1.007;
47+
mSamplingEnergyCorrectionPars[3] = 0.0667;
48+
mSamplingEnergyCorrectionPars[4] = -0.108;
49+
mSamplingEnergyCorrectionPars[5] = 0.0566;
50+
51+
mCrystalZCorrectionPars.reserve(9);
52+
mCrystalZCorrectionPars[0] = -0.005187;
53+
mCrystalZCorrectionPars[1] = 0.7301;
54+
mCrystalZCorrectionPars[2] = -0.7382;
55+
mCrystalZCorrectionPars[3] = 0.;
56+
mCrystalZCorrectionPars[4] = 0.;
57+
mCrystalZCorrectionPars[5] = 0.;
58+
mCrystalZCorrectionPars[6] = 0.;
59+
mCrystalZCorrectionPars[7] = 0.;
60+
mCrystalZCorrectionPars[8] = 0.;
61+
62+
mSamplingZCorrectionPars.reserve(9);
63+
mSamplingZCorrectionPars[0] = -2.137;
64+
mSamplingZCorrectionPars[1] = 6.400;
65+
mSamplingZCorrectionPars[2] = -3.342;
66+
mSamplingZCorrectionPars[3] = -0.1364;
67+
mSamplingZCorrectionPars[4] = 0.4019;
68+
mSamplingZCorrectionPars[5] = -0.1969;
69+
mSamplingZCorrectionPars[6] = 0.008223;
70+
mSamplingZCorrectionPars[7] = -0.02425;
71+
mSamplingZCorrectionPars[8] = 0.01190;
72+
3473
fCrystalShowerShape = new TF1("fCrystal", "x<[1] ? [0]*exp([3]*x+[4]*x*x+[5]*x*x*x) : (x<[2] ? [0]*[6]*exp([7]*x+[8]*x*x) : [0]*[9]*exp([10]*x+[11]*x*x))", 0, 15);
3574
double pc[12];
3675
pc[0] = 1. / 13.15;
@@ -354,25 +393,17 @@ void Clusterizer::evalClusters(std::vector<Cluster>& clusters)
354393

355394
// correct cluster energy and z position
356395
float eta = std::abs(cluster.getEta());
357-
float eCor = 1;
358-
float zCor = 0;
359396
bool isCrystal = geo.isCrystal(cluster.getDigitTowerId(0));
360-
if (isCrystal) {
361-
eCor = 0.00444 * std::pow(ee, -1.322) + (1.021 + 0.0018 * eta);
362-
if (mApplyCorrectionE)
363-
ee *= eCor;
364-
if (mApplyCorrectionZ)
365-
zCor = (-0.00518682 + 0.730052 * eta - 0.73817 * eta * eta);
366-
} else {
367-
eCor = 0.0033 * std::pow(ee, -2.09) + (1.007 + 0.0667 * eta - 0.108 * eta * eta + 0.0566 * eta * eta * eta);
368-
if (mApplyCorrectionE)
369-
ee *= eCor;
370-
if (mApplyCorrectionZ)
371-
zCor = (-2.13679 + 6.40009 * eta - 3.34233 * eta * eta) + (-0.136425 + 0.401887 * eta - 0.196851 * eta * eta) * ee + (0.00822276 - 0.0242512 * eta + 0.0118986 * eta * eta) * ee * ee;
397+
if (mApplyCorrectionE) {
398+
std::vector<double>& pe = isCrystal ? mCrystalEnergyCorrectionPars : mSamplingEnergyCorrectionPars;
399+
ee *= pe[0] * std::pow(ee, pe[1]) + pe[2] + pe[3] * eta + pe[4] * eta * eta + pe[5] * eta * eta * eta;
400+
cluster.setE(ee);
401+
}
402+
if (mApplyCorrectionZ) {
403+
std::vector<double>& pz = isCrystal ? mCrystalZCorrectionPars : mSamplingZCorrectionPars;
404+
float zCor = (pz[0] + pz[1] * eta + pz[2] * eta * eta) + (pz[3] + pz[4] * eta + pz[5] * eta * eta) * ee + (pz[6] + pz[7] * eta + pz[8] * eta * eta) * ee * ee;
405+
cluster.setZ(z > 0 ? z - zCor : z + zCor);
372406
}
373-
374-
cluster.setE(ee);
375-
cluster.setZ(cluster.getZ() - zCor);
376407

377408
// check if cluster is at the edge of detector module
378409
bool isEdge = 0;
@@ -385,7 +416,7 @@ void Clusterizer::evalClusters(std::vector<Cluster>& clusters)
385416
}
386417
cluster.setEdgeFlag(isEdge);
387418

388-
LOGF(debug, "Cluster coordinates: (%6.2f,%6.2f,%6.2f), eCor=%6.2f zCor=%6.2f", cluster.getX(), cluster.getY(), cluster.getZ(), eCor, zCor);
419+
LOGF(debug, "Cluster coordinates: (%6.2f,%6.2f,%6.2f)", cluster.getX(), cluster.getY(), cluster.getZ());
389420
}
390421
}
391422

0 commit comments

Comments
 (0)