File indexing completed on 2022-05-31 02:17:44
0001
0002
0003
0004
0005 template <typename Tf, typename Ti, typename TfLL1, typename Tf11, typename TfLLL>
0006 static inline void helixAtRFromIterativeCCS_impl(const Tf& __restrict__ inPar,
0007 const Ti& __restrict__ inChg,
0008 const Tf11& __restrict__ msRad,
0009 TfLL1& __restrict__ outPar,
0010 TfLLL& __restrict__ errorProp,
0011 Ti& __restrict__ outFailFlag,
0012 const int nmin,
0013 const int nmax,
0014 const int N_proc,
0015 const PropagationFlags pf) {
0016
0017
0018 #pragma omp simd
0019 for (int n = nmin; n < nmax; ++n) {
0020
0021 errorProp(n, 0, 0) = 1.f;
0022 errorProp(n, 1, 1) = 1.f;
0023 errorProp(n, 2, 2) = 1.f;
0024 errorProp(n, 3, 3) = 1.f;
0025 errorProp(n, 4, 4) = 1.f;
0026 errorProp(n, 5, 5) = 1.f;
0027 }
0028 float r0[nmax - nmin];
0029 #pragma omp simd
0030 for (int n = nmin; n < nmax; ++n) {
0031
0032 r0[n - nmin] = hipo(inPar(n, 0, 0), inPar(n, 1, 0));
0033 }
0034 float k[nmax - nmin];
0035 if (pf.use_param_b_field) {
0036 #pragma omp simd
0037 for (int n = nmin; n < nmax; ++n) {
0038 k[n - nmin] = inChg(n, 0, 0) * 100.f / (-Const::sol * Config::bFieldFromZR(inPar(n, 2, 0), r0[n - nmin]));
0039 }
0040 } else {
0041 #pragma omp simd
0042 for (int n = nmin; n < nmax; ++n) {
0043 k[n - nmin] = inChg(n, 0, 0) * 100.f / (-Const::sol * Config::Bfield);
0044 }
0045 }
0046 float r[nmax - nmin];
0047 #pragma omp simd
0048 for (int n = nmin; n < nmax; ++n) {
0049 r[n - nmin] = msRad(n, 0, 0);
0050 }
0051 float xin[nmax - nmin];
0052 float yin[nmax - nmin];
0053 float ipt[nmax - nmin];
0054 float phiin[nmax - nmin];
0055 float theta[nmax - nmin];
0056 #pragma omp simd
0057 for (int n = nmin; n < nmax; ++n) {
0058
0059
0060
0061
0062
0063 xin[n - nmin] = inPar(n, 0, 0);
0064 yin[n - nmin] = inPar(n, 1, 0);
0065 ipt[n - nmin] = inPar(n, 3, 0);
0066 phiin[n - nmin] = inPar(n, 4, 0);
0067 theta[n - nmin] = inPar(n, 5, 0);
0068
0069
0070 }
0071
0072 for (int n = nmin; n < nmax; ++n) {
0073 dprint_np(n,
0074 "input parameters"
0075 << " inPar(n, 0, 0)=" << std::setprecision(9) << inPar(n, 0, 0) << " inPar(n, 1, 0)="
0076 << std::setprecision(9) << inPar(n, 1, 0) << " inPar(n, 2, 0)=" << std::setprecision(9)
0077 << inPar(n, 2, 0) << " inPar(n, 3, 0)=" << std::setprecision(9) << inPar(n, 3, 0)
0078 << " inPar(n, 4, 0)=" << std::setprecision(9) << inPar(n, 4, 0)
0079 << " inPar(n, 5, 0)=" << std::setprecision(9) << inPar(n, 5, 0));
0080 }
0081
0082 float kinv[nmax - nmin];
0083 float pt[nmax - nmin];
0084 #pragma omp simd
0085 for (int n = nmin; n < nmax; ++n) {
0086 kinv[n - nmin] = 1.f / k[n - nmin];
0087 pt[n - nmin] = 1.f / ipt[n - nmin];
0088 }
0089 float D[nmax - nmin];
0090 float cosa[nmax - nmin];
0091 float sina[nmax - nmin];
0092 float cosah[nmax - nmin];
0093 float sinah[nmax - nmin];
0094 float id[nmax - nmin];
0095
0096 #pragma omp simd
0097 for (int n = nmin; n < nmax; ++n) {
0098 D[n - nmin] = 0.;
0099 }
0100
0101
0102 float cosPorT[nmax - nmin];
0103 float sinPorT[nmax - nmin];
0104 #pragma omp simd
0105 for (int n = nmin; n < nmax; ++n) {
0106 cosPorT[n - nmin] = std::cos(phiin[n - nmin]);
0107 }
0108 #pragma omp simd
0109 for (int n = nmin; n < nmax; ++n) {
0110 sinPorT[n - nmin] = std::sin(phiin[n - nmin]);
0111 }
0112
0113 float pxin[nmax - nmin];
0114 float pyin[nmax - nmin];
0115 #pragma omp simd
0116 for (int n = nmin; n < nmax; ++n) {
0117 pxin[n - nmin] = cosPorT[n - nmin] * pt[n - nmin];
0118 pyin[n - nmin] = sinPorT[n - nmin] * pt[n - nmin];
0119 }
0120
0121 for (int n = nmin; n < nmax; ++n) {
0122 dprint_np(n,
0123 "k=" << std::setprecision(9) << k[n - nmin] << " pxin=" << std::setprecision(9) << pxin[n - nmin]
0124 << " pyin=" << std::setprecision(9) << pyin[n - nmin] << " cosPorT=" << std::setprecision(9)
0125 << cosPorT[n - nmin] << " sinPorT=" << std::setprecision(9) << sinPorT[n - nmin]
0126 << " pt=" << std::setprecision(9) << pt[n - nmin]);
0127 }
0128
0129 float dDdx[nmax - nmin];
0130 float dDdy[nmax - nmin];
0131 float dDdipt[nmax - nmin];
0132 float dDdphi[nmax - nmin];
0133
0134 #pragma omp simd
0135 for (int n = nmin; n < nmax; ++n) {
0136 dDdipt[n - nmin] = 0.;
0137 dDdphi[n - nmin] = 0.;
0138 }
0139 #pragma omp simd
0140 for (int n = nmin; n < nmax; ++n) {
0141
0142 dDdx[n - nmin] = r0[n - nmin] > 0.f ? -xin[n - nmin] / r0[n - nmin] : 0.f;
0143 }
0144
0145 #pragma omp simd
0146 for (int n = nmin; n < nmax; ++n) {
0147 dDdy[n - nmin] = r0[n - nmin] > 0.f ? -yin[n - nmin] / r0[n - nmin] : 0.f;
0148 }
0149
0150 float oodotp[nmax - nmin];
0151 float x[nmax - nmin];
0152 float y[nmax - nmin];
0153 float oor0[nmax - nmin];
0154 float dadipt[nmax - nmin];
0155 float dadx[nmax - nmin];
0156 float dady[nmax - nmin];
0157 float pxca[nmax - nmin];
0158 float pxsa[nmax - nmin];
0159 float pyca[nmax - nmin];
0160 float pysa[nmax - nmin];
0161 float tmp[nmax - nmin];
0162 float tmpx[nmax - nmin];
0163 float tmpy[nmax - nmin];
0164 float pxinold[nmax - nmin];
0165
0166 CMS_UNROLL_LOOP_COUNT(Config::Niter)
0167 for (int i = 0; i < Config::Niter; ++i) {
0168 #pragma omp simd
0169 for (int n = nmin; n < nmax; ++n) {
0170
0171 r0[n - nmin] = hipo(outPar(n, 0, 0), outPar(n, 1, 0));
0172 }
0173
0174
0175
0176
0177
0178
0179
0180 #pragma omp simd
0181 for (int n = nmin; n < nmax; ++n) {
0182 oodotp[n - nmin] =
0183 r0[n - nmin] * pt[n - nmin] / (pxin[n - nmin] * outPar(n, 0, 0) + pyin[n - nmin] * outPar(n, 1, 0));
0184 }
0185
0186 #pragma omp simd
0187 for (int n = nmin; n < nmax; ++n) {
0188 if (oodotp[n - nmin] > 5.0f || oodotp[n - nmin] < 0)
0189 {
0190 outFailFlag(n, 0, 0) = 1;
0191 }
0192 }
0193
0194 #pragma omp simd
0195 for (int n = nmin; n < nmax; ++n) {
0196
0197
0198 id[n - nmin] =
0199 (oodotp[n - nmin] > 5.0f || oodotp[n - nmin] < 0) ? 0.0f : (r[n - nmin] - r0[n - nmin]) * oodotp[n - nmin];
0200 }
0201
0202 #pragma omp simd
0203 for (int n = nmin; n < nmax; ++n) {
0204 D[n - nmin] += id[n - nmin];
0205 }
0206
0207 if constexpr (Config::useTrigApprox) {
0208 #if !defined(__INTEL_COMPILER)
0209 #pragma omp simd
0210 #endif
0211 for (int n = nmin; n < nmax; ++n) {
0212 sincos4(id[n - nmin] * ipt[n - nmin] * kinv[n - nmin] * 0.5f, sinah[n - nmin], cosah[n - nmin]);
0213 }
0214 } else {
0215 #if !defined(__INTEL_COMPILER)
0216 #pragma omp simd
0217 #endif
0218 for (int n = nmin; n < nmax; ++n) {
0219 cosah[n - nmin] = std::cos(id[n - nmin] * ipt[n - nmin] * kinv[n - nmin] * 0.5f);
0220 sinah[n - nmin] = std::sin(id[n - nmin] * ipt[n - nmin] * kinv[n - nmin] * 0.5f);
0221 }
0222 }
0223
0224 #pragma omp simd
0225 for (int n = nmin; n < nmax; ++n) {
0226 cosa[n - nmin] = 1.f - 2.f * sinah[n - nmin] * sinah[n - nmin];
0227 sina[n - nmin] = 2.f * sinah[n - nmin] * cosah[n - nmin];
0228 }
0229
0230 for (int n = nmin; n < nmax; ++n) {
0231 dprint_np(n,
0232 "Attempt propagation from r="
0233 << r0[n - nmin] << " to r=" << r[n - nmin] << std::endl
0234 << " x=" << xin[n - nmin] << " y=" << yin[n - nmin] << " z=" << inPar(n, 2, 0)
0235 << " px=" << pxin[n - nmin] << " py=" << pyin[n - nmin]
0236 << " pz=" << pt[n - nmin] * std::tan(theta[n - nmin]) << " q=" << inChg(n, 0, 0) << std::endl
0237 << " r=" << std::setprecision(9) << r[n - nmin] << " r0=" << std::setprecision(9) << r0[n - nmin]
0238 << " id=" << std::setprecision(9) << id[n - nmin] << " dr=" << std::setprecision(9)
0239 << r[n - nmin] - r0[n - nmin] << " cosa=" << cosa[n - nmin] << " sina=" << sina[n - nmin]);
0240 }
0241
0242
0243 if (i + 1 != Config::Niter) {
0244 #pragma omp simd
0245 for (int n = nmin; n < nmax; ++n) {
0246 x[n - nmin] = outPar(n, 0, 0);
0247 y[n - nmin] = outPar(n, 1, 0);
0248 }
0249 #pragma omp simd
0250 for (int n = nmin; n < nmax; ++n) {
0251 oor0[n - nmin] =
0252 (r0[n - nmin] > 0.f && std::abs(r[n - nmin] - r0[n - nmin]) < 0.0001f) ? 1.f / r0[n - nmin] : 0.f;
0253 }
0254 #pragma omp simd
0255 for (int n = nmin; n < nmax; ++n) {
0256 dadipt[n - nmin] = id[n - nmin] * kinv[n - nmin];
0257 dadx[n - nmin] = -x[n - nmin] * ipt[n - nmin] * kinv[n - nmin] * oor0[n - nmin];
0258 dady[n - nmin] = -y[n - nmin] * ipt[n - nmin] * kinv[n - nmin] * oor0[n - nmin];
0259 pxca[n - nmin] = pxin[n - nmin] * cosa[n - nmin];
0260 pxsa[n - nmin] = pxin[n - nmin] * sina[n - nmin];
0261 pyca[n - nmin] = pyin[n - nmin] * cosa[n - nmin];
0262 pysa[n - nmin] = pyin[n - nmin] * sina[n - nmin];
0263 tmpx[n - nmin] = k[n - nmin] * dadx[n - nmin];
0264 }
0265
0266 #pragma omp simd
0267 for (int n = nmin; n < nmax; ++n) {
0268 dDdx[n - nmin] -= (x[n - nmin] * (1.f + tmpx[n - nmin] * (pxca[n - nmin] - pysa[n - nmin])) +
0269 y[n - nmin] * tmpx[n - nmin] * (pyca[n - nmin] + pxsa[n - nmin])) *
0270 oor0[n - nmin];
0271 }
0272
0273 #pragma omp simd
0274 for (int n = nmin; n < nmax; ++n) {
0275 tmpy[n - nmin] = k[n - nmin] * dady[n - nmin];
0276 }
0277 #pragma omp simd
0278 for (int n = nmin; n < nmax; ++n) {
0279 dDdy[n - nmin] -= (x[n - nmin] * tmpy[n - nmin] * (pxca[n - nmin] - pysa[n - nmin]) +
0280 y[n - nmin] * (1.f + tmpy[n - nmin] * (pyca[n - nmin] + pxsa[n - nmin]))) *
0281 oor0[n - nmin];
0282 }
0283 #pragma omp simd
0284 for (int n = nmin; n < nmax; ++n) {
0285
0286 tmp[n - nmin] = dadipt[n - nmin] * ipt[n - nmin];
0287 }
0288 #pragma omp simd
0289 for (int n = nmin; n < nmax; ++n) {
0290 dDdipt[n - nmin] -= k[n - nmin] *
0291 (x[n - nmin] * (pxca[n - nmin] * tmp[n - nmin] - pysa[n - nmin] * tmp[n - nmin] -
0292 pyca[n - nmin] - pxsa[n - nmin] + pyin[n - nmin]) +
0293 y[n - nmin] * (pyca[n - nmin] * tmp[n - nmin] + pxsa[n - nmin] * tmp[n - nmin] -
0294 pysa[n - nmin] + pxca[n - nmin] - pxin[n - nmin])) *
0295 pt[n - nmin] * oor0[n - nmin];
0296 }
0297 #pragma omp simd
0298 for (int n = nmin; n < nmax; ++n) {
0299 dDdphi[n - nmin] += k[n - nmin] *
0300 (x[n - nmin] * (pysa[n - nmin] - pxin[n - nmin] + pxca[n - nmin]) -
0301 y[n - nmin] * (pxsa[n - nmin] - pyin[n - nmin] + pyca[n - nmin])) *
0302 oor0[n - nmin];
0303 }
0304 }
0305
0306 #pragma omp simd
0307 for (int n = nmin; n < nmax; ++n) {
0308
0309 outPar(n, 0, 0) = outPar(n, 0, 0) + 2.f * k[n - nmin] * sinah[n - nmin] *
0310 (pxin[n - nmin] * cosah[n - nmin] - pyin[n - nmin] * sinah[n - nmin]);
0311 outPar(n, 1, 0) = outPar(n, 1, 0) + 2.f * k[n - nmin] * sinah[n - nmin] *
0312 (pyin[n - nmin] * cosah[n - nmin] + pxin[n - nmin] * sinah[n - nmin]);
0313 pxinold[n - nmin] = pxin[n - nmin];
0314 pxin[n - nmin] = pxin[n - nmin] * cosa[n - nmin] - pyin[n - nmin] * sina[n - nmin];
0315 pyin[n - nmin] = pyin[n - nmin] * cosa[n - nmin] + pxinold[n - nmin] * sina[n - nmin];
0316 }
0317 for (int n = nmin; n < nmax; ++n) {
0318 dprint_np(n,
0319 "outPar(n, 0, 0)=" << outPar(n, 0, 0) << " outPar(n, 1, 0)=" << outPar(n, 1, 0)
0320 << " pxin=" << pxin[n - nmin] << " pyin=" << pyin[n - nmin]);
0321 }
0322 }
0323
0324 float alpha[nmax - nmin];
0325 float dadphi[nmax - nmin];
0326
0327 #pragma omp simd
0328 for (int n = nmin; n < nmax; ++n) {
0329 alpha[n - nmin] = D[n - nmin] * ipt[n - nmin] * kinv[n - nmin];
0330 dadx[n - nmin] = dDdx[n - nmin] * ipt[n - nmin] * kinv[n - nmin];
0331 dady[n - nmin] = dDdy[n - nmin] * ipt[n - nmin] * kinv[n - nmin];
0332 dadipt[n - nmin] = (ipt[n - nmin] * dDdipt[n - nmin] + D[n - nmin]) * kinv[n - nmin];
0333 dadphi[n - nmin] = dDdphi[n - nmin] * ipt[n - nmin] * kinv[n - nmin];
0334 }
0335
0336 if constexpr (Config::useTrigApprox) {
0337 #pragma omp simd
0338 for (int n = nmin; n < nmax; ++n) {
0339 sincos4(alpha[n - nmin], sina[n - nmin], cosa[n - nmin]);
0340 }
0341 } else {
0342 #pragma omp simd
0343 for (int n = nmin; n < nmax; ++n) {
0344 cosa[n - nmin] = std::cos(alpha[n - nmin]);
0345 sina[n - nmin] = std::sin(alpha[n - nmin]);
0346 }
0347 }
0348 #pragma omp simd
0349 for (int n = nmin; n < nmax; ++n) {
0350 errorProp(n, 0, 0) = 1.f + k[n - nmin] * dadx[n - nmin] *
0351 (cosPorT[n - nmin] * cosa[n - nmin] - sinPorT[n - nmin] * sina[n - nmin]) *
0352 pt[n - nmin];
0353 errorProp(n, 0, 1) = k[n - nmin] * dady[n - nmin] *
0354 (cosPorT[n - nmin] * cosa[n - nmin] - sinPorT[n - nmin] * sina[n - nmin]) * pt[n - nmin];
0355 errorProp(n, 0, 2) = 0.f;
0356 errorProp(n, 0, 3) =
0357 k[n - nmin] *
0358 (cosPorT[n - nmin] * (ipt[n - nmin] * dadipt[n - nmin] * cosa[n - nmin] - sina[n - nmin]) +
0359 sinPorT[n - nmin] * ((1.f - cosa[n - nmin]) - ipt[n - nmin] * dadipt[n - nmin] * sina[n - nmin])) *
0360 pt[n - nmin] * pt[n - nmin];
0361 errorProp(n, 0, 4) =
0362 k[n - nmin] *
0363 (cosPorT[n - nmin] * dadphi[n - nmin] * cosa[n - nmin] - sinPorT[n - nmin] * dadphi[n - nmin] * sina[n - nmin] -
0364 sinPorT[n - nmin] * sina[n - nmin] + cosPorT[n - nmin] * cosa[n - nmin] - cosPorT[n - nmin]) *
0365 pt[n - nmin];
0366 errorProp(n, 0, 5) = 0.f;
0367 }
0368
0369 #pragma omp simd
0370 for (int n = nmin; n < nmax; ++n) {
0371 errorProp(n, 1, 0) = k[n - nmin] * dadx[n - nmin] *
0372 (sinPorT[n - nmin] * cosa[n - nmin] + cosPorT[n - nmin] * sina[n - nmin]) * pt[n - nmin];
0373 errorProp(n, 1, 1) = 1.f + k[n - nmin] * dady[n - nmin] *
0374 (sinPorT[n - nmin] * cosa[n - nmin] + cosPorT[n - nmin] * sina[n - nmin]) *
0375 pt[n - nmin];
0376 errorProp(n, 1, 2) = 0.f;
0377 errorProp(n, 1, 3) =
0378 k[n - nmin] *
0379 (sinPorT[n - nmin] * (ipt[n - nmin] * dadipt[n - nmin] * cosa[n - nmin] - sina[n - nmin]) +
0380 cosPorT[n - nmin] * (ipt[n - nmin] * dadipt[n - nmin] * sina[n - nmin] - (1.f - cosa[n - nmin]))) *
0381 pt[n - nmin] * pt[n - nmin];
0382 errorProp(n, 1, 4) =
0383 k[n - nmin] *
0384 (sinPorT[n - nmin] * dadphi[n - nmin] * cosa[n - nmin] + cosPorT[n - nmin] * dadphi[n - nmin] * sina[n - nmin] +
0385 sinPorT[n - nmin] * cosa[n - nmin] + cosPorT[n - nmin] * sina[n - nmin] - sinPorT[n - nmin]) *
0386 pt[n - nmin];
0387 errorProp(n, 1, 5) = 0.f;
0388 }
0389
0390 #pragma omp simd
0391 for (int n = nmin; n < nmax; ++n) {
0392
0393 cosPorT[n - nmin] = std::cos(theta[n - nmin]);
0394 }
0395
0396 #pragma omp simd
0397 for (int n = nmin; n < nmax; ++n) {
0398 sinPorT[n - nmin] = std::sin(theta[n - nmin]);
0399 }
0400
0401 #pragma omp simd
0402 for (int n = nmin; n < nmax; ++n) {
0403
0404 sinPorT[n - nmin] = 1.f / sinPorT[n - nmin];
0405 }
0406
0407 #pragma omp simd
0408 for (int n = nmin; n < nmax; ++n) {
0409 outPar(n, 2, 0) =
0410 inPar(n, 2, 0) + k[n - nmin] * alpha[n - nmin] * cosPorT[n - nmin] * pt[n - nmin] * sinPorT[n - nmin];
0411 errorProp(n, 2, 0) = k[n - nmin] * cosPorT[n - nmin] * dadx[n - nmin] * pt[n - nmin] * sinPorT[n - nmin];
0412 errorProp(n, 2, 1) = k[n - nmin] * cosPorT[n - nmin] * dady[n - nmin] * pt[n - nmin] * sinPorT[n - nmin];
0413 errorProp(n, 2, 2) = 1.f;
0414 errorProp(n, 2, 3) = k[n - nmin] * cosPorT[n - nmin] * (ipt[n - nmin] * dadipt[n - nmin] - alpha[n - nmin]) *
0415 pt[n - nmin] * pt[n - nmin] * sinPorT[n - nmin];
0416 errorProp(n, 2, 4) = k[n - nmin] * dadphi[n - nmin] * cosPorT[n - nmin] * pt[n - nmin] * sinPorT[n - nmin];
0417 errorProp(n, 2, 5) = -k[n - nmin] * alpha[n - nmin] * pt[n - nmin] * sinPorT[n - nmin] * sinPorT[n - nmin];
0418 }
0419
0420 #pragma omp simd
0421 for (int n = nmin; n < nmax; ++n) {
0422 outPar(n, 3, 0) = ipt[n - nmin];
0423 errorProp(n, 3, 0) = 0.f;
0424 errorProp(n, 3, 1) = 0.f;
0425 errorProp(n, 3, 2) = 0.f;
0426 errorProp(n, 3, 3) = 1.f;
0427 errorProp(n, 3, 4) = 0.f;
0428 errorProp(n, 3, 5) = 0.f;
0429 }
0430
0431 #pragma omp simd
0432 for (int n = nmin; n < nmax; ++n) {
0433 outPar(n, 4, 0) = inPar(n, 4, 0) + alpha[n - nmin];
0434 errorProp(n, 4, 0) = dadx[n - nmin];
0435 errorProp(n, 4, 1) = dady[n - nmin];
0436 errorProp(n, 4, 2) = 0.f;
0437 errorProp(n, 4, 3) = dadipt[n - nmin];
0438 errorProp(n, 4, 4) = 1.f + dadphi[n - nmin];
0439 errorProp(n, 4, 5) = 0.f;
0440 }
0441
0442 #pragma omp simd
0443 for (int n = nmin; n < nmax; ++n) {
0444 outPar(n, 5, 0) = theta[n - nmin];
0445 errorProp(n, 5, 0) = 0.f;
0446 errorProp(n, 5, 1) = 0.f;
0447 errorProp(n, 5, 2) = 0.f;
0448 errorProp(n, 5, 3) = 0.f;
0449 errorProp(n, 5, 4) = 0.f;
0450 errorProp(n, 5, 5) = 1.f;
0451 }
0452
0453 for (int n = nmin; n < nmax; ++n) {
0454 dprint_np(n,
0455 "propagation end, dump parameters"
0456 << std::endl
0457 << " pos = " << outPar(n, 0, 0) << " " << outPar(n, 1, 0) << " " << outPar(n, 2, 0) << "\t\t r="
0458 << std::sqrt(outPar(n, 0, 0) * outPar(n, 0, 0) + outPar(n, 1, 0) * outPar(n, 1, 0)) << std::endl
0459 << " mom = " << std::cos(outPar(n, 4, 0)) / outPar(n, 3, 0) << " "
0460 << std::sin(outPar(n, 4, 0)) / outPar(n, 3, 0) << " " << 1. / (outPar(n, 3, 0) * tan(outPar(n, 5, 0)))
0461 << "\t\tpT=" << 1. / std::abs(outPar(n, 3, 0)) << std::endl);
0462 }
0463
0464 #ifdef DEBUG
0465 for (int n = nmin; n < nmax; ++n) {
0466 if (n < N_proc) {
0467 dmutex_guard;
0468 std::cout << n << ": jacobian" << std::endl;
0469 printf("%5f %5f %5f %5f %5f %5f\n",
0470 errorProp(n, 0, 0),
0471 errorProp(n, 0, 1),
0472 errorProp(n, 0, 2),
0473 errorProp(n, 0, 3),
0474 errorProp(n, 0, 4),
0475 errorProp(n, 0, 5));
0476 printf("%5f %5f %5f %5f %5f %5f\n",
0477 errorProp(n, 1, 0),
0478 errorProp(n, 1, 1),
0479 errorProp(n, 1, 2),
0480 errorProp(n, 1, 3),
0481 errorProp(n, 1, 4),
0482 errorProp(n, 1, 5));
0483 printf("%5f %5f %5f %5f %5f %5f\n",
0484 errorProp(n, 2, 0),
0485 errorProp(n, 2, 1),
0486 errorProp(n, 2, 2),
0487 errorProp(n, 2, 3),
0488 errorProp(n, 2, 4),
0489 errorProp(n, 2, 5));
0490 printf("%5f %5f %5f %5f %5f %5f\n",
0491 errorProp(n, 3, 0),
0492 errorProp(n, 3, 1),
0493 errorProp(n, 3, 2),
0494 errorProp(n, 3, 3),
0495 errorProp(n, 3, 4),
0496 errorProp(n, 3, 5));
0497 printf("%5f %5f %5f %5f %5f %5f\n",
0498 errorProp(n, 4, 0),
0499 errorProp(n, 4, 1),
0500 errorProp(n, 4, 2),
0501 errorProp(n, 4, 3),
0502 errorProp(n, 4, 4),
0503 errorProp(n, 4, 5));
0504 printf("%5f %5f %5f %5f %5f %5f\n",
0505 errorProp(n, 5, 0),
0506 errorProp(n, 5, 1),
0507 errorProp(n, 5, 2),
0508 errorProp(n, 5, 3),
0509 errorProp(n, 5, 4),
0510 errorProp(n, 5, 5));
0511 printf("\n");
0512 }
0513 }
0514 #endif
0515 }