LArSoft  v09_90_00
Liquid Argon Software toolkit - https://larsoft.org/
LArRotationalTransformationPlugin.cc
Go to the documentation of this file.
1 
9 #include "Geometry/LArTPC.h"
10 
11 #include "Helpers/XmlHelper.h"
12 
13 #include "Managers/GeometryManager.h"
14 
15 #include "Objects/CartesianVector.h"
16 
17 #include "Pandora/Pandora.h"
18 
19 #include "larpandoracontent/LArPlugins/LArRotationalTransformationPlugin.h"
20 
21 #include <cmath>
22 
23 namespace lar_content
24 {
25 
26 using namespace pandora;
27 
28 LArRotationalTransformationPlugin::LArRotationalTransformationPlugin() :
29  m_thetaU(0.),
30  m_thetaV(0.),
31  m_thetaW(0.),
32  m_sinU(0.),
33  m_sinV(0.),
34  m_sinW(0.),
35  m_cosU(0.),
36  m_cosV(0.),
37  m_cosW(0.),
38  m_sinVminusU(0.),
39  m_sinWminusV(0.),
40  m_sinUminusW(0.),
41  m_maxAngularDiscrepancyU(0.03),
42  m_maxAngularDiscrepancyV(0.03),
43  m_maxAngularDiscrepancyW(0.03),
44  m_maxSigmaDiscrepancy(0.01)
45 {
46 }
47 
48 //------------------------------------------------------------------------------------------------------------------------------------------
49 
50 double LArRotationalTransformationPlugin::UVtoW(const double u, const double v) const
51 {
52  return (-1. * (u * m_sinWminusV + v * m_sinUminusW) / m_sinVminusU);
53 }
54 
55 //------------------------------------------------------------------------------------------------------------------------------------------
56 
57 double LArRotationalTransformationPlugin::VWtoU(const double v, const double w) const
58 {
59  return (-1. * (v * m_sinUminusW + w * m_sinVminusU) / m_sinWminusV);
60 }
61 
62 //------------------------------------------------------------------------------------------------------------------------------------------
63 
64 double LArRotationalTransformationPlugin::WUtoV(const double w, const double u) const
65 {
66  return (-1. * (u * m_sinWminusV + w * m_sinVminusU) / m_sinUminusW);
67 }
68 
69 //------------------------------------------------------------------------------------------------------------------------------------------
70 
71 double LArRotationalTransformationPlugin::UVtoY(const double u, const double v) const
72 {
73  return ((u * m_cosV - v * m_cosU) / m_sinVminusU);
74 }
75 
76 //------------------------------------------------------------------------------------------------------------------------------------------
77 
78 double LArRotationalTransformationPlugin::UVtoZ(const double u, const double v) const
79 {
80  return ((u * m_sinV - v * m_sinU) / m_sinVminusU);
81 }
82 
83 //------------------------------------------------------------------------------------------------------------------------------------------
84 
85 double LArRotationalTransformationPlugin::UWtoY(const double u, const double w) const
86 {
87  return ((w * m_cosU - u * m_cosW) / m_sinUminusW);
88 }
89 
90 //------------------------------------------------------------------------------------------------------------------------------------------
91 
92 double LArRotationalTransformationPlugin::UWtoZ(const double u, const double w) const
93 {
94  return ((w * m_sinU - u * m_sinW) / m_sinUminusW);
95 }
96 
97 //------------------------------------------------------------------------------------------------------------------------------------------
98 
99 double LArRotationalTransformationPlugin::VWtoY(const double v, const double w) const
100 {
101  return ((v * m_cosW - w * m_cosV) / m_sinWminusV);
102 }
103 
104 //------------------------------------------------------------------------------------------------------------------------------------------
105 
106 double LArRotationalTransformationPlugin::VWtoZ(const double v, const double w) const
107 {
108  return ((v * m_sinW - w * m_sinV) / m_sinWminusV);
109 }
110 
111 //------------------------------------------------------------------------------------------------------------------------------------------
112 
113 double LArRotationalTransformationPlugin::YZtoU(const double y, const double z) const
114 {
115  return (z * m_cosU - y * m_sinU);
116 }
117 
118 //------------------------------------------------------------------------------------------------------------------------------------------
119 
120 double LArRotationalTransformationPlugin::YZtoV(const double y, const double z) const
121 {
122  return (z * m_cosV - y * m_sinV);
123 }
124 
125 //------------------------------------------------------------------------------------------------------------------------------------------
126 
127 double LArRotationalTransformationPlugin::YZtoW(const double y, const double z) const
128 {
129  return (z * m_cosW - y * m_sinW);
130 }
131 
132 //------------------------------------------------------------------------------------------------------------------------------------------
133 
134 void LArRotationalTransformationPlugin::GetMinChiSquaredYZ(const double u, const double v, const double w, const double sigmaU,
135  const double sigmaV, const double sigmaW, double &y, double &z, double &chiSquared) const
136 {
137  const double sigmaU2(sigmaU * sigmaU), sigmaV2(sigmaV * sigmaV), sigmaW2(sigmaW * sigmaW);
138 
139  // Obtain expression for chi2, differentiate wrt y and z, set both results to zero and solve simultaneously. Here just paste-in result.
140  y = (sigmaW2 * v * m_cosU * m_cosV * m_sinU - sigmaW2 * u * m_cosV * m_cosV * m_sinU + sigmaV2 * w * m_cosU * m_cosW * m_sinU -
141  sigmaV2 * u * m_cosW * m_cosW * m_sinU - sigmaW2 * v * m_cosU * m_cosU * m_sinV + sigmaW2 * u * m_cosU * m_cosV * m_sinV +
142  sigmaU2 * w * m_cosV * m_cosW * m_sinV - sigmaU2 * v * m_cosW * m_cosW * m_sinV - sigmaV2 * w * m_cosU * m_cosU * m_sinW -
143  sigmaU2 * w * m_cosV * m_cosV * m_sinW + sigmaV2 * u * m_cosU * m_cosW * m_sinW + sigmaU2 * v * m_cosV * m_cosW * m_sinW) /
144  (sigmaW2 * m_cosV * m_cosV * m_sinU * m_sinU + sigmaV2 * m_cosW * m_cosW * m_sinU * m_sinU - 2. * sigmaW2 * m_cosU * m_cosV * m_sinU * m_sinV +
145  sigmaW2 * m_cosU * m_cosU * m_sinV * m_sinV + sigmaU2 * m_cosW * m_cosW * m_sinV * m_sinV -
146  2. * sigmaV2 * m_cosU * m_cosW * m_sinU * m_sinW - 2. * sigmaU2 * m_cosV * m_cosW * m_sinV * m_sinW +
147  sigmaV2 * m_cosU * m_cosU * m_sinW * m_sinW + sigmaU2 * m_cosV * m_cosV * m_sinW * m_sinW);
148 
149  z = (sigmaW2 * v * m_cosV * m_sinU * m_sinU + sigmaV2 * w * m_cosW * m_sinU * m_sinU - sigmaW2 * v * m_cosU * m_sinU * m_sinV -
150  sigmaW2 * u * m_cosV * m_sinU * m_sinV + sigmaW2 * u * m_cosU * m_sinV * m_sinV + sigmaU2 * w * m_cosW * m_sinV * m_sinV -
151  sigmaV2 * w * m_cosU * m_sinU * m_sinW - sigmaV2 * u * m_cosW * m_sinU * m_sinW - sigmaU2 * w * m_cosV * m_sinV * m_sinW -
152  sigmaU2 * v * m_cosW * m_sinV * m_sinW + sigmaV2 * u * m_cosU * m_sinW * m_sinW + sigmaU2 * v * m_cosV * m_sinW * m_sinW) /
153  (sigmaW2 * m_cosV * m_cosV * m_sinU * m_sinU + sigmaV2 * m_cosW * m_cosW * m_sinU * m_sinU - 2. * sigmaW2 * m_cosU * m_cosV * m_sinU * m_sinV +
154  sigmaW2 * m_cosU * m_cosU * m_sinV * m_sinV + sigmaU2 * m_cosW * m_cosW * m_sinV * m_sinV -
155  2. * sigmaV2 * m_cosU * m_cosW * m_sinU * m_sinW - 2. * sigmaU2 * m_cosV * m_cosW * m_sinV * m_sinW +
156  sigmaV2 * m_cosU * m_cosU * m_sinW * m_sinW + sigmaU2 * m_cosV * m_cosV * m_sinW * m_sinW);
157 
158  const double deltaU(u - LArRotationalTransformationPlugin::YZtoU(y, z));
159  const double deltaV(v - LArRotationalTransformationPlugin::YZtoV(y, z));
160  const double deltaW(w - LArRotationalTransformationPlugin::YZtoW(y, z));
161  chiSquared = ((deltaU * deltaU) / sigmaU2) + ((deltaV * deltaV) / sigmaV2) + ((deltaW * deltaW) / sigmaW2);
162 }
163 
164 //------------------------------------------------------------------------------------------------------------------------------------------
165 
166 void LArRotationalTransformationPlugin::GetMinChiSquaredYZ(const double u, const double v, const double w, const double sigmaU, const double sigmaV,
167  const double sigmaW, const double uFit, const double vFit, const double wFit, const double sigmaFit, double &y, double &z, double &chiSquared) const
168 {
169  const double sigmaU2(sigmaU * sigmaU), sigmaV2(sigmaV * sigmaV), sigmaW2(sigmaW * sigmaW), sigmaFit2(sigmaFit * sigmaFit);
170 
171  // Obtain expression for chi2, differentiate wrt y and z, set both results to zero and solve simultaneously. Here just paste-in result.
172  y = (vFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosV * m_sinU + vFit * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosV * m_sinU +
173  sigmaU2 * sigmaW2 * sigmaFit2 * v * m_cosU * m_cosV * m_sinU + sigmaW2 * sigmaFit2 * sigmaFit2 * v * m_cosU * m_cosV * m_sinU -
174  uFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosV * m_sinU - uFit * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosV * m_sinU -
175  sigmaV2 * sigmaW2 * sigmaFit2 * u * m_cosV * m_cosV * m_sinU - sigmaW2 * sigmaFit2 * sigmaFit2 * u * m_cosV * m_cosV * m_sinU +
176  wFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosW * m_sinU + wFit * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosW * m_sinU +
177  sigmaU2 * sigmaV2 * sigmaFit2 * w * m_cosU * m_cosW * m_sinU + sigmaV2 * sigmaFit2 * sigmaFit2 * w * m_cosU * m_cosW * m_sinU -
178  uFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_cosW * m_sinU - uFit * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosW * m_cosW * m_sinU -
179  sigmaV2 * sigmaW2 * sigmaFit2 * u * m_cosW * m_cosW * m_sinU - sigmaV2 * sigmaFit2 * sigmaFit2 * u * m_cosW * m_cosW * m_sinU -
180  vFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosU * m_sinV - vFit * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosU * m_sinV -
181  sigmaU2 * sigmaW2 * sigmaFit2 * v * m_cosU * m_cosU * m_sinV - sigmaW2 * sigmaFit2 * sigmaFit2 * v * m_cosU * m_cosU * m_sinV +
182  uFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosV * m_sinV + uFit * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosV * m_sinV +
183  sigmaV2 * sigmaW2 * sigmaFit2 * u * m_cosU * m_cosV * m_sinV + sigmaW2 * sigmaFit2 * sigmaFit2 * u * m_cosU * m_cosV * m_sinV +
184  wFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosW * m_sinV + wFit * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosW * m_sinV +
185  sigmaU2 * sigmaV2 * sigmaFit2 * w * m_cosV * m_cosW * m_sinV + sigmaU2 * sigmaFit2 * sigmaFit2 * w * m_cosV * m_cosW * m_sinV -
186  vFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_cosW * m_sinV - vFit * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosW * m_cosW * m_sinV -
187  sigmaU2 * sigmaW2 * sigmaFit2 * v * m_cosW * m_cosW * m_sinV - sigmaU2 * sigmaFit2 * sigmaFit2 * v * m_cosW * m_cosW * m_sinV -
188  wFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosU * m_sinW - wFit * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosU * m_sinW -
189  sigmaU2 * sigmaV2 * sigmaFit2 * w * m_cosU * m_cosU * m_sinW - sigmaV2 * sigmaFit2 * sigmaFit2 * w * m_cosU * m_cosU * m_sinW -
190  wFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosV * m_sinW - wFit * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosV * m_sinW -
191  sigmaU2 * sigmaV2 * sigmaFit2 * w * m_cosV * m_cosV * m_sinW - sigmaU2 * sigmaFit2 * sigmaFit2 * w * m_cosV * m_cosV * m_sinW +
192  uFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosW * m_sinW + uFit * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosU * m_cosW * m_sinW +
193  sigmaV2 * sigmaW2 * sigmaFit2 * u * m_cosU * m_cosW * m_sinW + sigmaV2 * sigmaFit2 * sigmaFit2 * u * m_cosU * m_cosW * m_sinW +
194  vFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosW * m_sinW + vFit * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosV * m_cosW * m_sinW +
195  sigmaU2 * sigmaW2 * sigmaFit2 * v * m_cosV * m_cosW * m_sinW + sigmaU2 * sigmaFit2 * sigmaFit2 * v * m_cosV * m_cosW * m_sinW) /
196  (sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosV * m_sinU * m_sinU + sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosV * m_sinU * m_sinU +
197  sigmaV2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosV * m_sinU * m_sinU + sigmaW2 * sigmaFit2 * sigmaFit2 * m_cosV * m_cosV * m_sinU * m_sinU +
198  sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_cosW * m_sinU * m_sinU + sigmaU2 * sigmaV2 * sigmaFit2 * m_cosW * m_cosW * m_sinU * m_sinU +
199  sigmaV2 * sigmaW2 * sigmaFit2 * m_cosW * m_cosW * m_sinU * m_sinU + sigmaV2 * sigmaFit2 * sigmaFit2 * m_cosW * m_cosW * m_sinU * m_sinU -
200  2. * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosV * m_sinU * m_sinV - 2. * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosV * m_sinU * m_sinV -
201  2. * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosV * m_sinU * m_sinV -
202  2. * sigmaW2 * sigmaFit2 * sigmaFit2 * m_cosU * m_cosV * m_sinU * m_sinV + sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosU * m_sinV * m_sinV +
203  sigmaU2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosU * m_sinV * m_sinV + sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosU * m_sinV * m_sinV +
204  sigmaW2 * sigmaFit2 * sigmaFit2 * m_cosU * m_cosU * m_sinV * m_sinV + sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_cosW * m_sinV * m_sinV +
205  sigmaU2 * sigmaV2 * sigmaFit2 * m_cosW * m_cosW * m_sinV * m_sinV + sigmaU2 * sigmaW2 * sigmaFit2 * m_cosW * m_cosW * m_sinV * m_sinV +
206  sigmaU2 * sigmaFit2 * sigmaFit2 * m_cosW * m_cosW * m_sinV * m_sinV - 2. * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosW * m_sinU * m_sinW -
207  2. * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosU * m_cosW * m_sinU * m_sinW - 2. * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosW * m_sinU * m_sinW -
208  2. * sigmaV2 * sigmaFit2 * sigmaFit2 * m_cosU * m_cosW * m_sinU * m_sinW -
209  2. * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosW * m_sinV * m_sinW - 2. * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosV * m_cosW * m_sinV * m_sinW -
210  2. * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosW * m_sinV * m_sinW -
211  2. * sigmaU2 * sigmaFit2 * sigmaFit2 * m_cosV * m_cosW * m_sinV * m_sinW +
212  sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosU * m_sinW * m_sinW + sigmaU2 * sigmaV2 * sigmaFit2 * m_cosU * m_cosU * m_sinW * m_sinW +
213  sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosU * m_sinW * m_sinW + sigmaV2 * sigmaFit2 * sigmaFit2 * m_cosU * m_cosU * m_sinW * m_sinW +
214  sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosV * m_sinW * m_sinW + sigmaU2 * sigmaV2 * sigmaFit2 * m_cosV * m_cosV * m_sinW * m_sinW +
215  sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosV * m_sinW * m_sinW + sigmaU2 * sigmaFit2 * sigmaFit2 * m_cosV * m_cosV * m_sinW * m_sinW);
216 
217  z = (vFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_sinU * m_sinU + vFit * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosV * m_sinU * m_sinU +
218  sigmaU2 * sigmaW2 * sigmaFit2 * v * m_cosV * m_sinU * m_sinU + sigmaW2 * sigmaFit2 * sigmaFit2 * v * m_cosV * m_sinU * m_sinU +
219  wFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_sinU * m_sinU + wFit * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosW * m_sinU * m_sinU +
220  sigmaU2 * sigmaV2 * sigmaFit2 * w * m_cosW * m_sinU * m_sinU + sigmaV2 * sigmaFit2 * sigmaFit2 * w * m_cosW * m_sinU * m_sinU -
221  vFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_sinU * m_sinV - vFit * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_sinU * m_sinV -
222  sigmaU2 * sigmaW2 * sigmaFit2 * v * m_cosU * m_sinU * m_sinV - sigmaW2 * sigmaFit2 * sigmaFit2 * v * m_cosU * m_sinU * m_sinV -
223  uFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_sinU * m_sinV - uFit * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_sinU * m_sinV -
224  sigmaV2 * sigmaW2 * sigmaFit2 * u * m_cosV * m_sinU * m_sinV - sigmaW2 * sigmaFit2 * sigmaFit2 * u * m_cosV * m_sinU * m_sinV +
225  uFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_sinV * m_sinV + uFit * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosU * m_sinV * m_sinV +
226  sigmaV2 * sigmaW2 * sigmaFit2 * u * m_cosU * m_sinV * m_sinV + sigmaW2 * sigmaFit2 * sigmaFit2 * u * m_cosU * m_sinV * m_sinV +
227  wFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_sinV * m_sinV + wFit * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosW * m_sinV * m_sinV +
228  sigmaU2 * sigmaV2 * sigmaFit2 * w * m_cosW * m_sinV * m_sinV + sigmaU2 * sigmaFit2 * sigmaFit2 * w * m_cosW * m_sinV * m_sinV -
229  wFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_sinU * m_sinW - wFit * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_sinU * m_sinW -
230  sigmaU2 * sigmaV2 * sigmaFit2 * w * m_cosU * m_sinU * m_sinW - sigmaV2 * sigmaFit2 * sigmaFit2 * w * m_cosU * m_sinU * m_sinW -
231  uFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_sinU * m_sinW - uFit * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosW * m_sinU * m_sinW -
232  sigmaV2 * sigmaW2 * sigmaFit2 * u * m_cosW * m_sinU * m_sinW - sigmaV2 * sigmaFit2 * sigmaFit2 * u * m_cosW * m_sinU * m_sinW -
233  wFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_sinV * m_sinW - wFit * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_sinV * m_sinW -
234  sigmaU2 * sigmaV2 * sigmaFit2 * w * m_cosV * m_sinV * m_sinW - sigmaU2 * sigmaFit2 * sigmaFit2 * w * m_cosV * m_sinV * m_sinW -
235  vFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_sinV * m_sinW - vFit * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosW * m_sinV * m_sinW -
236  sigmaU2 * sigmaW2 * sigmaFit2 * v * m_cosW * m_sinV * m_sinW - sigmaU2 * sigmaFit2 * sigmaFit2 * v * m_cosW * m_sinV * m_sinW +
237  uFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_sinW * m_sinW + uFit * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosU * m_sinW * m_sinW +
238  sigmaV2 * sigmaW2 * sigmaFit2 * u * m_cosU * m_sinW * m_sinW + sigmaV2 * sigmaFit2 * sigmaFit2 * u * m_cosU * m_sinW * m_sinW +
239  vFit * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_sinW * m_sinW + vFit * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosV * m_sinW * m_sinW +
240  sigmaU2 * sigmaW2 * sigmaFit2 * v * m_cosV * m_sinW * m_sinW + sigmaU2 * sigmaFit2 * sigmaFit2 * v * m_cosV * m_sinW * m_sinW) /
241  (sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosV * m_sinU * m_sinU + sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosV * m_sinU * m_sinU +
242  sigmaV2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosV * m_sinU * m_sinU + sigmaW2 * sigmaFit2 * sigmaFit2 * m_cosV * m_cosV * m_sinU * m_sinU +
243  sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_cosW * m_sinU * m_sinU + sigmaU2 * sigmaV2 * sigmaFit2 * m_cosW * m_cosW * m_sinU * m_sinU +
244  sigmaV2 * sigmaW2 * sigmaFit2 * m_cosW * m_cosW * m_sinU * m_sinU + sigmaV2 * sigmaFit2 * sigmaFit2 * m_cosW * m_cosW * m_sinU * m_sinU -
245  2. * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosV * m_sinU * m_sinV - 2. * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosV * m_sinU * m_sinV -
246  2. * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosV * m_sinU * m_sinV -
247  2. * sigmaW2 * sigmaFit2 * sigmaFit2 * m_cosU * m_cosV * m_sinU * m_sinV + sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosU * m_sinV * m_sinV +
248  sigmaU2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosU * m_sinV * m_sinV + sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosU * m_sinV * m_sinV +
249  sigmaW2 * sigmaFit2 * sigmaFit2 * m_cosU * m_cosU * m_sinV * m_sinV + sigmaU2 * sigmaV2 * sigmaW2 * m_cosW * m_cosW * m_sinV * m_sinV +
250  sigmaU2 * sigmaV2 * sigmaFit2 * m_cosW * m_cosW * m_sinV * m_sinV + sigmaU2 * sigmaW2 * sigmaFit2 * m_cosW * m_cosW * m_sinV * m_sinV +
251  sigmaU2 * sigmaFit2 * sigmaFit2 * m_cosW * m_cosW * m_sinV * m_sinV - 2. * sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosW * m_sinU * m_sinW -
252  2. * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosU * m_cosW * m_sinU * m_sinW - 2. * sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosW * m_sinU * m_sinW -
253  2. * sigmaV2 * sigmaFit2 * sigmaFit2 * m_cosU * m_cosW * m_sinU * m_sinW -
254  2. * sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosW * m_sinV * m_sinW - 2. * sigmaU2 * sigmaV2 * sigmaFit2 * m_cosV * m_cosW * m_sinV * m_sinW -
255  2. * sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosW * m_sinV * m_sinW -
256  2. * sigmaU2 * sigmaFit2 * sigmaFit2 * m_cosV * m_cosW * m_sinV * m_sinW +
257  sigmaU2 * sigmaV2 * sigmaW2 * m_cosU * m_cosU * m_sinW * m_sinW + sigmaU2 * sigmaV2 * sigmaFit2 * m_cosU * m_cosU * m_sinW * m_sinW +
258  sigmaV2 * sigmaW2 * sigmaFit2 * m_cosU * m_cosU * m_sinW * m_sinW + sigmaV2 * sigmaFit2 * sigmaFit2 * m_cosU * m_cosU * m_sinW * m_sinW +
259  sigmaU2 * sigmaV2 * sigmaW2 * m_cosV * m_cosV * m_sinW * m_sinW + sigmaU2 * sigmaV2 * sigmaFit2 * m_cosV * m_cosV * m_sinW * m_sinW +
260  sigmaU2 * sigmaW2 * sigmaFit2 * m_cosV * m_cosV * m_sinW * m_sinW + sigmaU2 * sigmaFit2 * sigmaFit2 * m_cosV * m_cosV * m_sinW * m_sinW);
261 
262  const double outputU(LArRotationalTransformationPlugin::YZtoU(y, z));
263  const double outputV(LArRotationalTransformationPlugin::YZtoV(y, z));
264  const double outputW(LArRotationalTransformationPlugin::YZtoW(y, z));
265 
266  const double deltaU(u - outputU), deltaV(v - outputV), deltaW(w - outputW);
267  const double deltaUFit(uFit - outputU), deltaVFit(vFit - outputV), deltaWFit(wFit - outputW);
268 
269  chiSquared = ((deltaU * deltaU) / sigmaU2) + ((deltaV * deltaV) / sigmaV2) + ((deltaW * deltaW) / sigmaW2) +
270  ((deltaUFit * deltaUFit) / sigmaFit2) + ((deltaVFit * deltaVFit) / sigmaFit2) + ((deltaWFit * deltaWFit) / sigmaFit2);
271 }
272 
273 //------------------------------------------------------------------------------------------------------------------------------------------
274 
275 StatusCode LArRotationalTransformationPlugin::Initialize()
276 {
277  const LArTPCMap &larTPCMap(this->GetPandora().GetGeometry()->GetLArTPCMap());
278 
279  if (larTPCMap.empty())
280  {
281  std::cout << "LArRotationalTransformationPlugin::Initialize - LArTPC description not registered with Pandora as required " << std::endl;
282  return STATUS_CODE_NOT_INITIALIZED;
283  }
284 
285  const LArTPC *const pFirstLArTPC(larTPCMap.begin()->second);
286  m_thetaU = pFirstLArTPC->GetWireAngleU();
287  m_thetaV = pFirstLArTPC->GetWireAngleV();
288  m_thetaW = pFirstLArTPC->GetWireAngleW();
289  const double sigmaUVW(pFirstLArTPC->GetSigmaUVW());
290 
291  m_sinU = std::sin(m_thetaU);
292  m_sinV = std::sin(m_thetaV);
293  m_sinW = std::sin(m_thetaW);
294  m_cosU = std::cos(m_thetaU);
295  m_cosV = std::cos(m_thetaV);
296  m_cosW = std::cos(m_thetaW);
297 
298  m_sinVminusU = std::sin(m_thetaV - m_thetaU);
299  m_sinWminusV = std::sin(m_thetaW - m_thetaV);
300  m_sinUminusW = std::sin(m_thetaU - m_thetaW);
301 
302  if ((std::fabs(m_sinVminusU) < std::numeric_limits<double>::epsilon()) ||
303  (std::fabs(m_sinWminusV) < std::numeric_limits<double>::epsilon()) || (std::fabs(m_sinUminusW) < std::numeric_limits<double>::epsilon()))
304  {
305  std::cout << "LArRotationalTransformationPlugin::Initialize - Equal wire angles; Plugin does not support provided LArTPC configurations. "
306  << std::endl;
307  return STATUS_CODE_INVALID_PARAMETER;
308  }
309 
310  for (const LArTPCMap::value_type &mapEntry : larTPCMap)
311  {
312  const LArTPC *const pLArTPC(mapEntry.second);
313 
314  if ((std::fabs(m_thetaU - pLArTPC->GetWireAngleU()) > m_maxAngularDiscrepancyU) ||
315  (std::fabs(m_thetaV - pLArTPC->GetWireAngleV()) > m_maxAngularDiscrepancyV) ||
316  (std::fabs(m_thetaW - pLArTPC->GetWireAngleW()) > m_maxAngularDiscrepancyW) ||
317  (std::fabs(sigmaUVW - pLArTPC->GetSigmaUVW()) > m_maxSigmaDiscrepancy))
318  {
319  std::cout << "LArRotationalTransformationPlugin::Initialize - Dissimilar drift volumes; Plugin does not support provided LArTPC configurations. "
320  << std::endl;
321  return STATUS_CODE_INVALID_PARAMETER;
322  }
323  }
324 
325  return STATUS_CODE_SUCCESS;
326 }
327 
328 //------------------------------------------------------------------------------------------------------------------------------------------
329 
330 pandora::StatusCode LArRotationalTransformationPlugin::ReadSettings(const pandora::TiXmlHandle xmlHandle)
331 {
332  PANDORA_RETURN_RESULT_IF_AND_IF(
333  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxAngularDiscrepancyU", m_maxAngularDiscrepancyU));
334 
335  PANDORA_RETURN_RESULT_IF_AND_IF(
336  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxAngularDiscrepancyV", m_maxAngularDiscrepancyV));
337 
338  PANDORA_RETURN_RESULT_IF_AND_IF(
339  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxAngularDiscrepancyW", m_maxAngularDiscrepancyW));
340 
341  PANDORA_RETURN_RESULT_IF_AND_IF(
342  STATUS_CODE_SUCCESS, STATUS_CODE_NOT_FOUND, !=, XmlHelper::ReadValue(xmlHandle, "MaxSigmaDiscrepancy", m_maxSigmaDiscrepancy));
343 
344  return STATUS_CODE_SUCCESS;
345 }
346 
347 } // namespace lar_content
lar_content::LArRotationalTransformationPlugin::GetMinChiSquaredYZ
virtual void GetMinChiSquaredYZ(const double u, const double v, const double w, const double sigmaU, const double sigmaV, const double sigmaW, double &y, double &z, double &chiSquared) const
Definition: LArRotationalTransformationPlugin.cc:134
lar_content::LArRotationalTransformationPlugin::m_cosV
double m_cosV
cos(thetaV)
Definition: LArRotationalTransformationPlugin.h:59
lar_content
Definition: CheatingBeamParticleIdTool.cc:18
lar_content::LArRotationalTransformationPlugin::UVtoY
virtual double UVtoY(const double u, const double v) const
Definition: LArRotationalTransformationPlugin.cc:71
lar_content::LArRotationalTransformationPlugin::m_sinU
double m_sinU
sin(thetaU)
Definition: LArRotationalTransformationPlugin.h:55
lar_content::LArRotationalTransformationPlugin::UVtoW
virtual double UVtoW(const double u, const double v) const
Definition: LArRotationalTransformationPlugin.cc:50
lar_content::LArRotationalTransformationPlugin::Initialize
pandora::StatusCode Initialize()
Definition: LArRotationalTransformationPlugin.cc:275
lar_content::LArRotationalTransformationPlugin::YZtoW
virtual double YZtoW(const double y, const double z) const
Definition: LArRotationalTransformationPlugin.cc:127
lar_content::LArRotationalTransformationPlugin::UVtoZ
virtual double UVtoZ(const double u, const double v) const
Definition: LArRotationalTransformationPlugin.cc:78
lar_content::LArRotationalTransformationPlugin::m_sinVminusU
double m_sinVminusU
sin(thetaV - thetaU)
Definition: LArRotationalTransformationPlugin.h:61
lar_content::LArRotationalTransformationPlugin::VWtoU
virtual double VWtoU(const double v, const double w) const
Definition: LArRotationalTransformationPlugin.cc:57
y
Float_t y
Definition: compare.C:6
z
Double_t z
Definition: plot.C:276
lar_content::LArRotationalTransformationPlugin::m_cosU
double m_cosU
cos(thetaU)
Definition: LArRotationalTransformationPlugin.h:58
lar_content::LArRotationalTransformationPlugin::m_thetaW
double m_thetaW
inclination of W wires (radians)
Definition: LArRotationalTransformationPlugin.h:53
lar_content::LArRotationalTransformationPlugin::m_maxAngularDiscrepancyU
double m_maxAngularDiscrepancyU
Maximum allowed difference between u wire angles between LArTPCs.
Definition: LArRotationalTransformationPlugin.h:65
lar_content::LArRotationalTransformationPlugin::m_sinUminusW
double m_sinUminusW
sin(thetaU - thetaW)
Definition: LArRotationalTransformationPlugin.h:63
lar_content::LArRotationalTransformationPlugin::UWtoZ
virtual double UWtoZ(const double u, const double w) const
Definition: LArRotationalTransformationPlugin.cc:92
lar_content::LArRotationalTransformationPlugin::m_sinW
double m_sinW
sin(thetaW)
Definition: LArRotationalTransformationPlugin.h:57
lar_content::LArRotationalTransformationPlugin::m_thetaV
double m_thetaV
inclination of V wires (radians)
Definition: LArRotationalTransformationPlugin.h:52
lar_content::LArRotationalTransformationPlugin::m_cosW
double m_cosW
cos(thetaW)
Definition: LArRotationalTransformationPlugin.h:60
LArRotationalTransformationPlugin.h
Header file for the rotational transformation plugin class.
lar_content::LArRotationalTransformationPlugin::YZtoV
virtual double YZtoV(const double y, const double z) const
Definition: LArRotationalTransformationPlugin.cc:120
lar_content::LArRotationalTransformationPlugin::m_sinV
double m_sinV
sin(thetaV)
Definition: LArRotationalTransformationPlugin.h:56
lar_content::LArRotationalTransformationPlugin::m_maxAngularDiscrepancyV
double m_maxAngularDiscrepancyV
Maximum allowed difference between v wire angles between LArTPCs.
Definition: LArRotationalTransformationPlugin.h:66
lar_content::LArRotationalTransformationPlugin::m_maxAngularDiscrepancyW
double m_maxAngularDiscrepancyW
Maximum allowed difference between w wire angles between LArTPCs.
Definition: LArRotationalTransformationPlugin.h:67
lar_content::LArRotationalTransformationPlugin::m_sinWminusV
double m_sinWminusV
sin(thetaW - thetaV)
Definition: LArRotationalTransformationPlugin.h:62
lar_content::LArRotationalTransformationPlugin::ReadSettings
pandora::StatusCode ReadSettings(const pandora::TiXmlHandle xmlHandle)
Definition: LArRotationalTransformationPlugin.cc:330
lar_content::LArRotationalTransformationPlugin::VWtoZ
virtual double VWtoZ(const double v, const double w) const
Definition: LArRotationalTransformationPlugin.cc:106
lar_content::LArRotationalTransformationPlugin::WUtoV
virtual double WUtoV(const double w, const double u) const
Definition: LArRotationalTransformationPlugin.cc:64
lar_content::LArRotationalTransformationPlugin::m_thetaU
double m_thetaU
inclination of U wires (radians)
Definition: LArRotationalTransformationPlugin.h:51
lar_content::LArRotationalTransformationPlugin::YZtoU
virtual double YZtoU(const double y, const double z) const
Definition: LArRotationalTransformationPlugin.cc:113
lar_content::LArRotationalTransformationPlugin::m_maxSigmaDiscrepancy
double m_maxSigmaDiscrepancy
Maximum allowed difference between like wire sigma values between LArTPCs.
Definition: LArRotationalTransformationPlugin.h:68
w
Float_t w
Definition: plot.C:20
lar_content::LArRotationalTransformationPlugin::LArRotationalTransformationPlugin
LArRotationalTransformationPlugin()
Default constructor.
Definition: LArRotationalTransformationPlugin.cc:28
lar_content::LArRotationalTransformationPlugin::VWtoY
virtual double VWtoY(const double v, const double w) const
Definition: LArRotationalTransformationPlugin.cc:99
lar_content::LArRotationalTransformationPlugin::UWtoY
virtual double UWtoY(const double u, const double w) const
Definition: LArRotationalTransformationPlugin.cc:85