CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
Functions
VVIObjDetails Namespace Reference

Functions

double cosint (double x)
 Private version of the cosine and sine integral. More...
 
template<typename F >
int dzero (double a, double b, double &x0, double &rv, double eps, int mxf, F func)
 
double expint (double x)
 Private version of the sine integral. More...
 
double f1 (double x, double const *h_)
 Private version of the exponential integral. More...
 
double f2 (double x, double const *h_)
 
void sincosint (double x, double &sint, double &cint)
 
double sinint (double x)
 Private version of the cosine integral. More...
 

Function Documentation

double VVIObjDetails::cosint ( double  x)

Private version of the cosine and sine integral.

Definition at line 195 of file VVIObj.cc.

References funct::cos(), alignCSCRings::e, create_public_lumi_plots::log, AlCaHLTBitMon_ParallelJobs::p, createTree::pp, lumiQueryAPI::q, funct::sin(), vdt::x, detailsBasic3DVector::y, and zero.

195  {
196  // Initialized data
197 
198  const double zero = 0.;
199  const double one = 1.;
200  const double two = 2.;
201  const double eight = 8.;
202  const double ce = .57721566490153;
203  const double c__[14] = { 1.9405491464836,.9413409132865,
204  -.579845034293,.3091572011159,-.0916101792208,.0164437407515,
205  -.0019713091952,1.692538851e-4,-1.09393296e-5,5.522386e-7,
206  -2.23995e-8,7.465e-10,-2.08e-11,5e-13 };
207  const double p[23] = { .96074783975204,-.0371138962124,
208  .00194143988899,-1.7165988425e-4,2.112637753e-5,-3.27163257e-6,
209  6.0069212e-7,-1.2586794e-7,2.932563e-8,-7.45696e-9,2.04105e-9,
210  -5.9502e-10,1.8323e-10,-5.921e-11,1.997e-11,-7e-12,2.54e-12,
211  -9.5e-13,3.7e-13,-1.4e-13,6e-14,-2e-14,1e-14 };
212  const double q[20] = { .98604065696238,-.0134717382083,
213  4.5329284117e-4,-3.067288652e-5,3.13199198e-6,-4.2110196e-7,
214  6.907245e-8,-1.318321e-8,2.83697e-9,-6.7329e-10,1.734e-10,
215  -4.787e-11,1.403e-11,-4.33e-12,1.4e-12,-4.7e-13,1.7e-13,-6e-14,
216  2e-14,-1e-14 };
217 
218  // System generated locals
219  double d__1;
220 
221  // Local variables
222  double h__;
223  int i__;
224  double r__, y, b0, b1, b2, pp, qq, alfa;
225 
226  // If x==0, return same
227 
228  if (x == zero) {
229  return zero;
230  }
231  if (fabs(x) <= eight) {
232  y = x / eight;
233  // Computing 2nd power
234  d__1 = y;
235  h__ = two * (d__1 * d__1) - one;
236  alfa = -two * h__;
237  b1 = zero;
238  b2 = zero;
239  for (i__ = 13; i__ >= 0; --i__) {
240  b0 = c__[i__] - alfa * b1 - b2;
241  b2 = b1;
242  b1 = b0;
243  }
244  b1 = ce + log((fabs(x))) - b0 + h__ * b2;
245  } else {
246  r__ = one / x;
247  y = eight * r__;
248  // Computing 2nd power
249  d__1 = y;
250  h__ = two * (d__1 * d__1) - one;
251  alfa = -two * h__;
252  b1 = zero;
253  b2 = zero;
254  for (i__ = 22; i__ >= 0; --i__) {
255  b0 = p[i__] - alfa * b1 - b2;
256  b2 = b1;
257  b1 = b0;
258  }
259  pp = b0 - h__ * b2;
260  b1 = zero;
261  b2 = zero;
262  for (i__ = 19; i__ >= 0; --i__) {
263  b0 = q[i__] - alfa * b1 - b2;
264  b2 = b1;
265  b1 = b0;
266  }
267  qq = b0 - h__ * b2;
268  b1 = r__ * (qq * sin(x) - r__ * pp * cos(x));
269  }
270  return b1;
271  } // cosint
tuple pp
Definition: createTree.py:15
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
x
Definition: VDTMath.h:216
template<typename F >
int VVIObjDetails::dzero ( double  a,
double  b,
double &  x0,
double &  rv,
double  eps,
int  mxf,
F  func 
)

Definition at line 572 of file VVIObj.cc.

References f, f1(), f2(), python.connectstrParser::f3, benchmark_cfg::fa, benchmark_cfg::fb, alignCSCRings::ff, max(), and min.

Referenced by VVIObj::VVIObj().

573  {
574  /* System generated locals */
575  double d__1, d__2, d__3, d__4;
576 
577  // Local variables
578  double f1, f2, f3, u1, u2, x1, x2, u3, u4, x3, ca, cb, cc, fa, fb, ee, ff;
579  int mc;
580  double xa, xb, fx, xx, su4;
581 
582  xa = std::min(a,b);
583  xb = std::max(a,b);
584  fa = func(xa);
585  fb = func(xb);
586  if (fa * fb > 0.) {
587  rv = (xb - xa) * -2;
588  x0 = 0.;
589  return 1;
590  }
591  mc = 0;
592  L1:
593  x0 = (xa + xb) * .5;
594  rv = x0 - xa;
595  ee = eps * (fabs(x0) + 1);
596  if (rv <= ee) {
597  rv = ee;
598  ff = func(x0);
599  return 0;
600  }
601  f1 = fa;
602  x1 = xa;
603  f2 = fb;
604  x2 = xb;
605  L2:
606  fx = func(x0);
607  ++mc;
608  if (mc > mxf) {
609  rv = (d__1 = xb - xa, fabs(d__1)) * -.5;
610  x0 = 0.;
611  return 0;
612  }
613  if (fx * fa > 0.) {
614  xa = x0;
615  fa = fx;
616  } else {
617  xb = x0;
618  fb = fx;
619  }
620  L3:
621  u1 = f1 - f2;
622  u2 = x1 - x2;
623  u3 = f2 - fx;
624  u4 = x2 - x0;
625  if (u2 == 0. || u4 == 0.) {goto L1;}
626  f3 = fx;
627  x3 = x0;
628  u1 /= u2;
629  u2 = u3 / u4;
630  ca = u1 - u2;
631  cb = (x1 + x2) * u2 - (x2 + x0) * u1;
632  cc = (x1 - x0) * f1 - x1 * (ca * x1 + cb);
633  if (ca == 0.) {
634  if (cb == 0.) {goto L1;}
635  x0 = -cc / cb;
636  } else {
637  u3 = cb / (ca * 2);
638  u4 = u3 * u3 - cc / ca;
639  if (u4 < 0.) {goto L1;}
640  su4 = fabs(u4);
641  if (x0 + u3 < 0.f) {su4 = -su4;}
642  x0 = -u3 + su4;
643  }
644  if (x0 < xa || x0 > xb) {goto L1;}
645  // Computing MIN
646  d__3 = (d__1 = x0 - x3, fabs(d__1)), d__4 = (d__2 = x0 - x2, fabs(d__2));
647  rv = std::min(d__3,d__4);
648  ee = eps * (fabs(x0) + 1);
649  if (rv > ee) {
650  f1 = f2;
651  x1 = x2;
652  f2 = f3;
653  x2 = x3;
654  goto L2;
655  }
656  fx = func(x0);
657  if (fx == 0.) {
658  rv = ee;
659  ff = func(x0);
660  return 0;
661  }
662  if (fx * fa < 0.) {
663  xx = x0 - ee;
664  if (xx <= xa) {
665  rv = ee;
666  ff = func(x0);
667  return 0;
668  }
669  ff = func(xx);
670  fb = ff;
671  xb = xx;
672  } else {
673  xx = x0 + ee;
674  if (xx >= xb) {
675  rv = ee;
676  ff = func(x0);
677  return 0;
678  }
679  ff = func(xx);
680  fa = ff;
681  xa = xx;
682  }
683  if (fx * ff > 0.) {
684  mc += 2;
685  if (mc > mxf) {
686  rv = (d__1 = xb - xa, fabs(d__1)) * -.5;
687  x0 = 0.;
688  return 0;
689  }
690  f1 = f3;
691  x1 = x3;
692  f2 = fx;
693  x2 = x0;
694  x0 = xx;
695  fx = ff;
696  goto L3;
697  }
698  /* L4: */
699  rv = ee;
700  ff = func(x0);
701  return 0;
702  } // dzero
#define min(a, b)
Definition: mlp_lapack.h:161
const T & max(const T &a, const T &b)
double f[11][100]
double b
Definition: hdecay.h:120
double a
Definition: hdecay.h:121
double VVIObjDetails::expint ( double  x)

Private version of the sine integral.

Definition at line 450 of file VVIObj.cc.

References create_public_lumi_plots::exp, create_public_lumi_plots::log, p1, p2, p3, p4, q1, q2, v, vdt::x, detailsBasic3DVector::y, and zero.

Referenced by f2(), and VVIObj::VVIObj().

450  {
451 
452  // Initialized data
453 
454  const double zero = 0.;
455  const double q2[7] = { .10340013040487,3.319092135933,
456  20.449478501379,41.280784189142,32.426421069514,10.041164382905,
457  1. };
458  const double p3[6] = { -2.3909964453136,-147.98219500504,
459  -254.3763397689,-119.55761038372,-19.630408535939,-.9999999999036
460  };
461  const double q3[6] = { 177.60070940351,530.68509610812,
462  462.23027156148,156.81843364539,21.630408494238,1. };
463  const double p4[8] = { -8.6693733995107,-549.14226552109,
464  -4210.0161535707,-249301.39345865,-119623.66934925,
465  -22174462.775885,3892804.213112,-391546073.8091 };
466  const double q4[8] = { 34.171875,-1607.0892658722,35730.029805851,
467  -483547.43616216,4285596.2461175,-24903337.574054,89192576.757561,
468  -165254299.72521 };
469  const double a1[8] = { -2.1808638152072,-21.901023385488,
470  9.3081638566217,25.076281129356,-33.184253199722,60.121799083008,
471  -43.253113287813,1.0044310922808 };
472  const double b1[8] = { 0.,3.9370770185272,300.89264837292,
473  -6.2504116167188,1003.6743951673,14.325673812194,2736.2411988933,
474  .52746885196291 };
475  const double a2[8] = { -3.4833465360285,-18.65454548834,
476  -8.2856199414064,-32.34673303054,17.960168876925,1.7565631546961,
477  -1.9502232128966,.99999429607471 };
478  const double b2[8] = { 0.,69.500065588743,57.283719383732,
479  25.777638423844,760.76114800773,28.951672792514,-3.4394226689987,
480  1.0008386740264 };
481  const double a3[6] = { -27.780928934438,-10.10479081576,
482  -9.1483008216736,-5.0223317461851,-3.0000077799358,
483  1.0000000000704 };
484  const double one = 1.;
485  const double b3[6] = { 0.,122.39993926823,2.7276100778779,
486  -7.1897518395045,-2.9990118065262,1.999999942826 };
487  const double two = 2.;
488  const double three = 3.;
489  const double x0 = .37250741078137;
490  const double xl[6] = { -24.,-12.,-6.,0.,1.,4. };
491  const double p1[5] = { 4.293125234321,39.894153870321,
492  292.52518866921,425.69682638592,-434.98143832952 };
493  const double q1[5] = { 1.,18.899288395003,150.95038744251,
494  568.05252718987,753.58564359843 };
495  const double p2[7] = { .43096783946939,6.9052252278444,
496  23.019255939133,24.378408879132,9.0416155694633,.99997957705159,
497  4.656271079751e-7 };
498 
499  /* Local variables */
500  double v, y, ap, bp, aq, dp, bq, dq;
501 
502  if (x <= xl[0]) {
503  ap = a3[0] - x;
504  for ( int i__ = 2; i__ <= 5; ++i__) {
505  /* L1: */
506  ap = a3[i__ - 1] - x + b3[i__ - 1] / ap;
507  }
508  y = exp(-x) / x * (one - (a3[5] + b3[5] / ap) / x);
509  } else if (x <= xl[1]) {
510  ap = a2[0] - x;
511  for ( int i__ = 2; i__ <= 7; ++i__) {
512  ap = a2[i__ - 1] - x + b2[i__ - 1] / ap;
513  }
514  y = exp(-x) / x * (a2[7] + b2[7] / ap);
515  } else if (x <= xl[2]) {
516  ap = a1[0] - x;
517  for ( int i__ = 2; i__ <= 7; ++i__) {
518  ap = a1[i__ - 1] - x + b1[i__ - 1] / ap;
519  }
520  y = exp(-x) / x * (a1[7] + b1[7] / ap);
521  } else if (x < xl[3]) {
522  v = -two * (x / three + one);
523  bp = zero;
524  dp = p4[0];
525  for ( int i__ = 2; i__ <= 8; ++i__) {
526  ap = bp;
527  bp = dp;
528  dp = p4[i__ - 1] - ap + v * bp;
529  }
530  bq = zero;
531  dq = q4[0];
532  for ( int i__ = 2; i__ <= 8; ++i__) {
533  aq = bq;
534  bq = dq;
535  dq = q4[i__ - 1] - aq + v * bq;
536  }
537  y = -log(-x / x0) + (x + x0) * (dp - ap) / (dq - aq);
538  } else if (x == xl[3]) {
539  return zero;
540  } else if (x < xl[4]) {
541  ap = p1[0];
542  aq = q1[0];
543  for ( int i__ = 2; i__ <= 5; ++i__) {
544  ap = p1[i__ - 1] + x * ap;
545  aq = q1[i__ - 1] + x * aq;
546  }
547  y = -log(x) + ap / aq;
548  } else if (x <= xl[5]) {
549  y = one / x;
550  ap = p2[0];
551  aq = q2[0];
552  for ( int i__ = 2; i__ <= 7; ++i__) {
553  ap = p2[i__ - 1] + y * ap;
554  aq = q2[i__ - 1] + y * aq;
555  }
556  y = exp(-x) * ap / aq;
557  } else {
558  y = one / x;
559  ap = p3[0];
560  aq = q3[0];
561  for ( int i__ = 2; i__ <= 6; ++i__) {
562  ap = p3[i__ - 1] + y * ap;
563  aq = q3[i__ - 1] + y * aq;
564  }
565  y = exp(-x) * y * (one + y * ap / aq);
566  }
567  return y;
568 } // expint
double q2[4]
Definition: TauolaWrapper.h:88
double p4[4]
Definition: TauolaWrapper.h:92
double p2[4]
Definition: TauolaWrapper.h:90
double q1[4]
Definition: TauolaWrapper.h:87
double p1[4]
Definition: TauolaWrapper.h:89
x
Definition: VDTMath.h:216
mathSSE::Vec4< T > v
double p3[4]
Definition: TauolaWrapper.h:91
double VVIObjDetails::f1 ( double  x,
double const *  h_ 
)
inline

Private version of the exponential integral.

Definition at line 33 of file VVIObj.cc.

References create_public_lumi_plots::log, and vdt::x.

Referenced by dzero(), and VVIObj::VVIObj().

33 { return h_[0]+h_[1]*std::log(h_[2]*x)-h_[3]*x;}
x
Definition: VDTMath.h:216
double VVIObjDetails::f2 ( double  x,
double const *  h_ 
)
inline

Definition at line 34 of file VVIObj.cc.

References abs, create_public_lumi_plots::exp, expint(), and create_public_lumi_plots::log.

Referenced by dzero(), and VVIObj::VVIObj().

34 { return h_[4]-x+h_[5]*(std::log(std::abs(x))+expint(x))-h_[6]*std::exp(-x);}
#define abs(x)
Definition: mlp_lapack.h:159
x
Definition: VDTMath.h:216
double expint(double x)
Private version of the sine integral.
Definition: VVIObj.cc:450
void VVIObjDetails::sincosint ( double  x,
double &  sint,
double &  cint 
)

Definition at line 347 of file VVIObj.cc.

References funct::cos(), alignCSCRings::e, create_public_lumi_plots::log, AlCaHLTBitMon_ParallelJobs::p, createTree::pp, lumiQueryAPI::q, funct::sin(), vdt::x, detailsBasic3DVector::y, and zero.

Referenced by VVIObj::VVIObj().

347  {
348  // Initialized data
349 
350  const double zero = 0.;
351  const double one = 1.;
352  const double two = 2.;
353  const double eight = 8.;
354  const double ce = .57721566490153;
355  const double pih = 1.5707963267949;
356  const double s__[14] = { 1.9522209759531,-.6884042321257,
357  .4551855132256,-.1804571236838,.0410422133759,-.0059586169556,
358  6.001427414e-4,-4.44708329e-5,2.5300782e-6,-1.141308e-7,4.1858e-9,
359  -1.273e-10,3.3e-12,-1e-13 };
360 
361  const double c__[14] = { 1.9405491464836,.9413409132865,
362  -.579845034293,.3091572011159,-.0916101792208,.0164437407515,
363  -.0019713091952,1.692538851e-4,-1.09393296e-5,5.522386e-7,
364  -2.23995e-8,7.465e-10,-2.08e-11,5e-13 };
365 
366  const double p[23] = { .96074783975204,-.0371138962124,
367  .00194143988899,-1.7165988425e-4,2.112637753e-5,-3.27163257e-6,
368  6.0069212e-7,-1.2586794e-7,2.932563e-8,-7.45696e-9,2.04105e-9,
369  -5.9502e-10,1.8323e-10,-5.921e-11,1.997e-11,-7e-12,2.54e-12,
370  -9.5e-13,3.7e-13,-1.4e-13,6e-14,-2e-14,1e-14 };
371  const double q[20] = { .98604065696238,-.0134717382083,
372  4.5329284117e-4,-3.067288652e-5,3.13199198e-6,-4.2110196e-7,
373  6.907245e-8,-1.318321e-8,2.83697e-9,-6.7329e-10,1.734e-10,
374  -4.787e-11,1.403e-11,-4.33e-12,1.4e-12,-4.7e-13,1.7e-13,-6e-14,
375  2e-14,-1e-14 };
376 
377  // System generated locals
378  double d__1;
379 
380  // Local variables
381  double h__;
382  int i__;
383  double r__, y, b0, b1, b2, pp, qq, alfa;
384 
385  sint=0;
386  cint=0;
387 
388 
389  if (fabs(x) <= eight) {
390  y = x / eight;
391  // Computing 2nd power
392  d__1 = y;
393  h__ = two * (d__1 * d__1) - one;
394  alfa = -two * h__;
395 
396  // cos
397  if (x!=0) {
398  b1 = zero;
399  b2 = zero;
400  for (i__ = 13; i__ >= 0; --i__) {
401  b0 = c__[i__] - alfa * b1 - b2;
402  b2 = b1;
403  b1 = b0;
404  }
405  cint = ce + log((fabs(x))) - b0 + h__ * b2;
406  }
407  // sin
408  b1 = zero;
409  b2 = zero;
410  for (i__ = 13; i__ >= 0; --i__) {
411  b0 = s__[i__] - alfa * b1 - b2;
412  b2 = b1;
413  b1 = b0;
414  }
415  sint = y * (b0 - b2);
416 
417  } else {
418  r__ = one / x;
419  y = eight * r__;
420  // Computing 2nd power
421  d__1 = y;
422  h__ = two * (d__1 * d__1) - one;
423  alfa = -two * h__;
424  b1 = zero;
425  b2 = zero;
426  for (i__ = 22; i__ >= 0; --i__) {
427  b0 = p[i__] - alfa * b1 - b2;
428  b2 = b1;
429  b1 = b0;
430  }
431  pp = b0 - h__ * b2;
432  b1 = zero;
433  b2 = zero;
434  for (i__ = 19; i__ >= 0; --i__) {
435  b0 = q[i__] - alfa * b1 - b2;
436  b2 = b1;
437  b1 = b0;
438  }
439  qq = b0 - h__ * b2;
440  // cos
441  cint = r__ * (qq * sin(x) - r__ * pp * cos(x));
442  // sin
443  d__1 = pih;
444  if(x < 0.) d__1 = -d__1;
445  sint = d__1 - r__ * (r__ * pp * sin(x) + qq * cos(x));
446  }
447  }
tuple pp
Definition: createTree.py:15
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
x
Definition: VDTMath.h:216
double VVIObjDetails::sinint ( double  x)

Private version of the cosine integral.

Definition at line 273 of file VVIObj.cc.

References funct::cos(), alignCSCRings::e, AlCaHLTBitMon_ParallelJobs::p, createTree::pp, lumiQueryAPI::q, alignCSCRings::s, funct::sin(), vdt::x, detailsBasic3DVector::y, and zero.

273  {
274  // Initialized data
275 
276  const double zero = 0.;
277  const double one = 1.;
278  const double two = 2.;
279  const double eight = 8.;
280  const double pih = 1.5707963267949;
281  const double s[14] = { 1.9522209759531,-.6884042321257,
282  .4551855132256,-.1804571236838,.0410422133759,-.0059586169556,
283  6.001427414e-4,-4.44708329e-5,2.5300782e-6,-1.141308e-7,4.1858e-9,
284  -1.273e-10,3.3e-12,-1e-13 };
285  const double p[23] = { .96074783975204,-.0371138962124,
286  .00194143988899,-1.7165988425e-4,2.112637753e-5,-3.27163257e-6,
287  6.0069212e-7,-1.2586794e-7,2.932563e-8,-7.45696e-9,2.04105e-9,
288  -5.9502e-10,1.8323e-10,-5.921e-11,1.997e-11,-7e-12,2.54e-12,
289  -9.5e-13,3.7e-13,-1.4e-13,6e-14,-2e-14,1e-14 };
290  const double q[20] = { .98604065696238,-.0134717382083,
291  4.5329284117e-4,-3.067288652e-5,3.13199198e-6,-4.2110196e-7,
292  6.907245e-8,-1.318321e-8,2.83697e-9,-6.7329e-10,1.734e-10,
293  -4.787e-11,1.403e-11,-4.33e-12,1.4e-12,-4.7e-13,1.7e-13,-6e-14,
294  2e-14,-1e-14 };
295 
296  // System generated locals
297  double d__1;
298 
299  // Local variables
300  double h__;
301  int i__;
302  double r__, y, b0, b1, b2, pp, qq, alfa;
303 
304  if (fabs(x) <= eight) {
305  y = x / eight;
306  d__1 = y;
307  h__ = two * (d__1 * d__1) - one;
308  alfa = -two * h__;
309  b1 = zero;
310  b2 = zero;
311  for (i__ = 13; i__ >= 0; --i__) {
312  b0 = s[i__] - alfa * b1 - b2;
313  b2 = b1;
314  b1 = b0;
315  }
316  b1 = y * (b0 - b2);
317  } else {
318  r__ = one / x;
319  y = eight * r__;
320  d__1 = y;
321  h__ = two * (d__1 * d__1) - one;
322  alfa = -two * h__;
323  b1 = zero;
324  b2 = zero;
325  for (i__ = 22; i__ >= 0; --i__) {
326  b0 = p[i__] - alfa * b1 - b2;
327  b2 = b1;
328  b1 = b0;
329  }
330  pp = b0 - h__ * b2;
331  b1 = zero;
332  b2 = zero;
333  for (i__ = 19; i__ >= 0; --i__) {
334  b0 = q[i__] - alfa * b1 - b2;
335  b2 = b1;
336  b1 = b0;
337  }
338  qq = b0 - h__ * b2;
339  d__1 = fabs(pih);
340  if(x < 0.) d__1 = -d__1;
341  b1 = d__1 - r__ * (r__ * pp * sin(x) + qq * cos(x));
342  }
343 
344  return b1;
345  } // sinint
tuple pp
Definition: createTree.py:15
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
x
Definition: VDTMath.h:216