31 std::map<unsigned int, AlignmentResult> &
results,
55 map<unsigned int, unsigned int> offset_map;
57 offsets.push_back(dim);
66 unsigned int detId = dit.first;
77 double r_xx, r_xy, r_xz;
78 double r_yx, r_yy, r_yz;
79 double r_zx, r_zy, r_zz;
80 rotation.GetComponents(r_xx, r_xy, r_xz, r_yx, r_yy, r_yz, r_zx, r_zy, r_zz);
83 throw cms::Exception(
"PPS") <<
"IdealResult::Solve: only rotations about z are supported.";
85 double rot_z = atan2(r_yx, r_xx);
98 const auto &
d =
geom.getDirectionData(1);
103 const auto &
d =
geom.getDirectionData(2);
110 chi_tr(offsets[qci] + idx) =
v;
115 vector<TVectorD> inaccessibleModes;
118 TVectorD fm_ShX_gl(dim);
120 TVectorD fm_ShX_lp(dim);
122 TVectorD fm_ShY_gl(dim);
124 TVectorD fm_ShY_lp(dim);
129 const auto &
geom = sp.second;
134 const double d2x =
geom.getDirectionData(2).dx;
135 const double d2y =
geom.getDirectionData(2).dy;
138 fm_ShX_gl(offset2 + qIndex) = d2x;
139 fm_ShX_lp(offset2 + qIndex) = d2x *
geom.z;
140 fm_ShY_gl(offset2 + qIndex) = d2y;
141 fm_ShY_lp(offset2 + qIndex) = d2y *
geom.z;
148 const double d1x =
geom.getDirectionData(1).dx;
149 const double d1y =
geom.getDirectionData(1).dy;
150 const double d2x =
geom.getDirectionData(2).dx;
151 const double d2y =
geom.getDirectionData(2).dy;
154 fm_ShX_gl(offset1 + qIndex1) = d1x;
155 fm_ShX_lp(offset1 + qIndex1) = d1x *
geom.z;
156 fm_ShY_gl(offset1 + qIndex1) = d1y;
157 fm_ShY_lp(offset1 + qIndex1) = d1y *
geom.z;
160 fm_ShX_gl(offset2 + qIndex2) = d2x;
161 fm_ShX_lp(offset2 + qIndex2) = d2x *
geom.z;
162 fm_ShY_gl(offset2 + qIndex2) = d2y;
163 fm_ShY_lp(offset2 + qIndex2) = d2y *
geom.z;
167 inaccessibleModes.push_back(fm_ShX_gl);
168 inaccessibleModes.push_back(fm_ShX_lp);
169 inaccessibleModes.push_back(fm_ShY_gl);
170 inaccessibleModes.push_back(fm_ShY_lp);
174 TVectorD fm_RotZ_gl(dim);
176 TVectorD fm_RotZ_lp(dim);
181 const auto &
geom = sp.second;
183 for (
int m = 0;
m < 2; ++
m) {
185 TVectorD *fmp =
nullptr;
200 const double as_x = -rho *
geom.sy;
201 const double as_y = +rho *
geom.sx;
204 const double d2x =
geom.getDirectionData(2).dx;
205 const double d2y =
geom.getDirectionData(2).dy;
209 fm(offset2 + qIndex2) = d2x * as_x + d2y * as_y;
213 const double d1x =
geom.getDirectionData(1).dx;
214 const double d1y =
geom.getDirectionData(1).dy;
215 const double d2x =
geom.getDirectionData(2).dx;
216 const double d2y =
geom.getDirectionData(2).dy;
220 fm(offset1 + qIndex1) = d1x * as_x + d1y * as_y;
224 fm(offset2 + qIndex2) = d2x * as_x + d2y * as_y;
229 inaccessibleModes.push_back(fm_RotZ_gl);
230 inaccessibleModes.push_back(fm_RotZ_lp);
234 TMatrixD
C(dim, constraints.size());
235 TVectorD
V(constraints.size());
236 for (
unsigned int i = 0;
i < constraints.size();
i++) {
237 V(
i) = constraints[
i].val;
241 const TVectorD &
cv = constraints[
i].coef.find(quantityClass)->second;
242 for (
int k = 0;
k < cv.GetNrows();
k++) {
249 TMatrixD
I(dim, inaccessibleModes.size());
250 for (
unsigned int i = 0;
i < inaccessibleModes.size(); ++
i) {
251 for (
int j = 0;
j < inaccessibleModes[
i].GetNrows(); ++
j)
252 I(
j,
i) = inaccessibleModes[
i](
j);
256 TMatrixD CT(TMatrixD::kTransposed,
C);
257 TMatrixD CTI(CT *
I);
259 const TMatrixD &
A = CTI;
260 TMatrixD AT(TMatrixD::kTransposed, A);
261 TMatrixD ATA(AT * A);
262 TMatrixD ATA_inv(TMatrixD::kInverted, ATA);
264 TVectorD
b =
V - CT * chi_tr;
266 TVectorD La(ATA_inv * AT * b);
268 TVectorD chi_exp(chi_tr + I * La);
280 const auto &
v = chi_exp(offsets[qci] + idx);
298 results[dit.first] =
r;
AlignmentTask * task
the tasked to be completed
void setShR2(const double &v)
const Translation & translation() const
void setShR1(const double &v)
IdealResult()
dummy constructor (not to be used)
std::vector< QuantityClass > quantityClasses
list of quantity classes to be optimized
Abstract parent for all (track-based) alignment algorithms.
unsigned int solve(const std::vector< AlignmentConstraint > &, std::map< unsigned int, AlignmentResult > &result, TDirectory *dir) override
void setShZ(const double &v)
bool resolveRotZ
whether to resolve detector rotations around z
Represents an alignment task.
uint32_t T const *__restrict__ uint32_t const *__restrict__ int32_t int Histo::index_type cudaStream_t V
detector shifts in first readout direction
void setRotZ(const double &v)
edm::ESHandle< CTPPSGeometry > gMisaligned
AlignmentGeometry geometry
the geometry for this task
unsigned int quantitiesOfClass(QuantityClass) const
returns the number of quantities of the given class
void begin(const CTPPSGeometry *geometryReal, const CTPPSGeometry *geometryMisaligned) override
prepare for processing
constexpr int subdetId() const
get the contents of the subdetector field (not cast into any detector's numbering enum) ...
Abs< T >::type abs(const T &t)
const std::complex< double > I
The manager class for TOTEM RP geometry.
detector shifts in second readout direction
const std::map< unsigned int, DetGeometry > & getSensorMap() const
edm::ESHandle< CTPPSGeometry > gReal
bool resolveShR
whether to resolve detector shifts in readout direction(s)
detector rotations around z
Result of CTPPS track-based alignment.
const DetGeometry & get(unsigned int id) const
retrieves sensor geometry
Base class for CTPPS detector IDs.
signed int getQuantityIndex(QuantityClass cl, unsigned int detId) const
returns measurement index (if non-existent, returns -1)
static unsigned int const shift
const RotationMatrix & rotation() const