Constructor.
Constructor Set Vavilov parameters kappa and beta2 and define whether to calculate density fcn or distribution fcn
- Parameters:
-
| kappa | - (input) Vavilov kappa parameter [0.01 (Landau-like) < kappa < 10. (Gaussian-like)] |
| beta2 | - (input) Vavilov beta2 parameter (square of particle speed in v/c units) |
| mode | - (input) set to 0 to calculate the density function and to 1 to calculate the distribution function |
Definition at line 48 of file VVIObj.cc.
References a_, abs, b_, alignmentValidation::c1, funct::cos(), VVIObjDetails::dzero(), alignCSCRings::e, create_public_lumi_plots::exp, VVIObjDetails::expint(), python::connectstrParser::f1, python::connectstrParser::f2, gen::k, prof2calltree::l, create_public_lumi_plots::log, mode_, n, omega_, lumiQueryAPI::q, q2, funct::sin(), VVIObjDetails::sincosint(), t0_, t1_, t_, x, and x0_.
: mode_(mode) {
const double xp[9] = { 9.29,2.47,.89,.36,.15,.07,.03,.02,0.0 };
const double xq[7] = { .012,.03,.08,.26,.87,3.83,11.0 };
double h_[7];
double q, u, x, c1, c2, c3, c4, d1, h4, h5, h6, q2, x1, d, ll, ul, xf1, xf2, rv;
int lp, lq, k, l, n;
if(kappa < 0.01) kappa = 0.01;
if(kappa > 10.) kappa = 10.;
if(beta2 < 0.) beta2 = 0.;
if(beta2 > 1.) beta2 = 1.;
h_[4] = 1. - beta2*0.42278433999999998 + 7.6/kappa;
h_[5] = beta2;
h_[6] = 1. - beta2;
h4 = -7.6/kappa - (beta2 * .57721566 + 1);
h5 = log(kappa);
h6 = 1./kappa;
t0_ = (h4 - h_[4]*h5 - (h_[4] + beta2)*(log(h_[4]) + VVIObjDetails::expint(h_[4])) + exp(-h_[4]))/h_[4];
for (lp = 0; lp < 9; ++lp) {
if (kappa >= xp[lp]) break;
}
ll = -lp - 1.5;
for (lq = 0; lq < 7; ++lq) {
if (kappa <= xq[lq]) break;
}
ul = lq - 6.5;
auto f2 = [h_](double x) { return h_[4]-x+h_[5]*(std::log(std::abs(x))+VVIObjDetails::expint(x))-h_[6]*std::exp(-x);};
VVIObjDetails::dzero(ll, ul, u, rv, 1.e-5, 1000, f2);
q = 1./u;
t1_ = h4 * q - h5 - (beta2 * q + 1) * (log((fabs(u))) + VVIObjDetails::expint(u)) + exp(-u) * q;
t_ = t1_ - t0_;
omega_ = 6.2831853000000004/t_;
h_[0] = kappa * (beta2 * .57721566 + 2.) + 9.9166128600000008;
if (kappa >= .07) {h_[0] += 6.90775527;}
h_[1] = beta2 * kappa;
h_[2] = h6 * omega_;
h_[3] = omega_ * 1.5707963250000001;
auto f1 = [h_](double x){ return h_[0]+h_[1]*std::log(h_[2]*x)-h_[3]*x;};
VVIObjDetails::dzero(5., 155., x0_, rv, 1.e-5, 1000, f1);
n = x0_ + 1.;
d = exp(kappa * (beta2 * (.57721566 - h5) + 1.)) * .31830988654751274;
a_[n - 1] = 0.;
if (mode_ == 0) {
a_[n - 1] = omega_ * .31830988654751274;
}
q = -1.;
q2 = 2.;
for (k = 1; k < n; ++k) {
l = n - k;
x = omega_ * k;
x1 = h6 * x;
VVIObjDetails::sincosint(x1,c2,c1);
c1 = log(x) - c1;
c3 = sin(x1);
c4 = cos(x1);
xf1 = kappa * (beta2 * c1 - c4) - x * c2;
xf2 = x * c1 + kappa * (c3 + beta2 * c2) + t0_ * x;
if (mode_ == 0) {
d1 = q * d * omega_ * exp(xf1);
a_[l - 1] = d1 * cos(xf2);
b_[l - 1] = -d1 * sin(xf2);
} else {
d1 = q * d * exp(xf1)/k;
a_[l - 1] = d1 * sin(xf2);
b_[l - 1] = d1 * cos(xf2);
a_[n - 1] += q2 * a_[l - 1];
}
q = -q;
q2 = -q2;
}
}
1.7.3