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 194 of file VVIObj.cc.

References funct::cos(), ExpressReco_HICollisions_FallBack::e, funct::log(), L1TEmulatorMonitor_cff::p, createTree::pp, lumiQueryAPI::q, funct::sin(), ExpressReco_HICollisions_FallBack::x, ExpressReco_HICollisions_FallBack::y, and zero.

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

Definition at line 571 of file VVIObj.cc.

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

Referenced by VVIObj::VVIObj().

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

References funct::exp(), funct::log(), p1, p2, p3, p4, q1, q2, v, ExpressReco_HICollisions_FallBack::x, ExpressReco_HICollisions_FallBack::y, and zero.

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

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

References funct::log(), and ExpressReco_HICollisions_FallBack::x.

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

32 { return h_[0]+h_[1]*std::log(h_[2]*x)-h_[3]*x;}
Log< T >::type log(const T &t)
Definition: Log.h:22
double VVIObjDetails::f2 ( double  x,
double const *  h_ 
)
inline

Definition at line 33 of file VVIObj.cc.

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

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

33 { 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
Exp< T >::type exp(const T &t)
Definition: Exp.h:22
Log< T >::type log(const T &t)
Definition: Log.h:22
double expint(double x)
Private version of the sine integral.
Definition: VVIObj.cc:449
void VVIObjDetails::sincosint ( double  x,
double &  sint,
double &  cint 
)

Definition at line 346 of file VVIObj.cc.

References funct::cos(), ExpressReco_HICollisions_FallBack::e, funct::log(), L1TEmulatorMonitor_cff::p, createTree::pp, lumiQueryAPI::q, funct::sin(), ExpressReco_HICollisions_FallBack::x, ExpressReco_HICollisions_FallBack::y, and zero.

Referenced by VVIObj::VVIObj().

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

Private version of the cosine integral.

Definition at line 272 of file VVIObj.cc.

References funct::cos(), ExpressReco_HICollisions_FallBack::e, L1TEmulatorMonitor_cff::p, createTree::pp, lumiQueryAPI::q, asciidump::s, funct::sin(), ExpressReco_HICollisions_FallBack::x, ExpressReco_HICollisions_FallBack::y, and zero.

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