/* cteq3.f -- translated by f2c (version 19961001). You must link the resulting object file with the libraries: -lf2c -lm (in that order) */ #include "f2c.h" /* Table of constant values */ static doublereal c_b20 = 1.122; static doublereal c_b21 = .9906; static doublereal c_b22 = 1.; static doublereal c_b33 = .9838; static doublereal c_b34 = 1.001; static doublereal c_b44 = 1.105; static doublereal c_b45 = .9818; static doublereal c_b46 = 1.002; /* Version 3 CTEQ distribution function in a parametrized form. */ /* By: H.L. Lai, J. Botts, J. Huston, J.G. Morfin, J.F. Owens, J. Qiu, */ /* W.K. Tung & H. Weerts; Preprint MSU-HEP/41024, CTEQ 404 */ /* This file contains three versions of the same CTEQ3 parton distributions: */ /* Two "front-end" subprograms: */ /* FUNCTION Ctq3Pd (Iset, Iparton, X, Q, Irt) */ /* returns the PROBABILITY density for a GIVEN flavor; */ /* SUBROUTINE Ctq3Pds(Iset, Pdf, XX, QQ, Irt) */ /* returns an array of MOMENTUM densities for ALL flavors; */ /* One lower-level subprogram: */ /* FUNCTION Ctq3df (Iset, Iprtn, XX, QQ, Irt) */ /* returns the MOMENTUM density of a GIVEN valence or sea distribution. */ /* One supplementary function to return the QCD lambda parameter */ /* concerning these distributions is also included (see below). */ /* Although DOUBLE PRECISION is used, conversion to SINGLE PRECISION */ /* is straightforward by removing the */ /* Implicit Double Precision statements. */ /* Since this is an initial distribution of version 3, it is */ /* useful for the authors to maintain a record of the distribution */ /* list in case there are revisions or corrections. */ /* In the interest of maintaining the integrity of this package, */ /* please do not freely distribute this program package; instead, refer*/ /* any interested colleagues to direct their request for a copy to: */ /* Lai@cteq11.pa.msu.edu or Tung@msupa.pa.msu.edu. */ /* If you have detailed questions concerning these CTEQ3 distributions, */ /* or if you find problems/bugs using this initial distribution, direct */ /* inquires to Hung-Liang Lai or Wu-Ki Tung. */ /* ------------------------------------------- */ /* Detailed instructions follow. */ /* Name convention for CTEQ distributions: CTEQnSx where */ /* n : version number (currently n = 3) */ /* S : factorization scheme label: = [M L D] for [MS-bar LO DIS] */ /* resp. */ /* x : special characteristics, if any */ /* (e.g. S(F) for singular (flat) small-x, L for "LEP lambda value")*/ /* (not applicable to CTEQ3 since only three standard sets are given.)*/ /* Explanation of functional arguments: */ /* Iset is the set label; in this version, Iset = 1, 2, 3 */ /* correspond to the following CTEQ global fits: */ /* cteq3M : best fit in the MS-bar scheme */ /* cteq3L : best fit in Leading order QCD */ /* cteq3D : best fit in the DIS scheme */ /* Iprtn is the parton label (6, 5, 4, 3, 2, 1, 0, -1, ......, -6) */ /* for (t, b, c, s, d, u, g, u_bar, ..., t_bar) */ /* *** WARNING: We use the parton label 2 as D-quark, and 1 as U-quark which*/ /* might be different with your labels. */ /* X, Q are the usual x, Q; */ /* Irt is a return error code (see individual modules for explanation). */ /* --------------------------------------------- */ /* Since the QCD Lambda value for the various sets are needed more often than */ /* the other parameters in most applications, a special function */ /* Wlamd3 (Iset, Iorder, Neff) is provided */ /* which returns the lambda value for Neff = 4,5,6 effective flavors as well a s*/ /* the order these values pertain to. */ /* ---------------------------------------------- */ /* The range of (x, Q) used in this round of global analysis is, approxi-*/ /* mately, 0.01 < x < 0.75 ; and 4 GeV^2 < Q^2 < 400 GeV^2 for fixed targe t*/ /* experiments and 0.0001 < x < 0.1 from HERA data. */ /* The range of (x, Q) used in the reparametrization of the QCD evolved */ /* parton distributions is 10E-6 < x < 1 ; 1.6 GeV < Q < 10 TeV. The */ /* functional form of this parametrization is: */ /* A0 * x^A1 * (1-x)^A2 * (1 + A3 * x^A4) * [log(1+1/x)]^A5 */ /* with the A'coefficients being smooth functions of Q. For heavy quarks,*/ /* a threshold factor is applied to A0 which simulates the proper Q-dependenc e*/ /* of the QCD evolution in that region according to the renormalization */ /* scheme defined in Collins-Tung, Nucl. Phys. B278, 934 (1986). */ /* Since this function is positive definite and smooth, it provides sensible */ /* extrapolations of the parton distributions if they are called beyond */ /* the original range in an application. There is no artificial boundaries*/ /* or sharp cutoff's. */ /* ------------------------------------------------ */ doublereal ctq3pd_(iset, iparton, x, q, irt) integer *iset, *iparton; doublereal *x, *q; integer *irt; { /* System generated locals */ integer i__1; doublereal ret_val; /* Local variables */ extern doublereal ctq3df_(); static doublereal vl, sea; static integer ifl, jfl; /* This function returns the CTEQ parton distributions f^Iset_Iprtn/proto n*/ /* --- the PROBABILITY density */ /* (Iset, Iparton, X, Q): explained above; */ /* Irt : return error code: see module Ctq3df for explanation. */ ifl = *iparton; jfl = abs(ifl); /* Valence */ if (ifl == 1 || ifl == 2) { vl = ctq3df_(iset, &ifl, x, q, irt); } else { vl = (float)0.; } /* Sea */ i__1 = -jfl; sea = ctq3df_(iset, &i__1, x, q, irt); /* Full (probability) Distribu tion*/ ret_val = (vl + sea) / *x; return ret_val; /* ************************* */ } /* ctq3pd_ */ /* Subroutine */ int ctq3pds_(iset, pdf, x, q, irt) integer *iset; doublereal *pdf, *x, *q; integer *irt; { static integer i__; extern doublereal ctq3df_(); static integer irt1; /* This function returns the CTEQ parton distributions xf^Iset_Iprtn/prot on*/ /* --- the Momentum density in array form */ /* (Iset, X, Q): explained in header comment lines; */ /* Irt : return error code -- cumulated over flavors: */ /* see module Ctq3df for explanation on individual flavors. */ /* Pdf (Iparton); */ /* Iparton = -6, -5, ...0, 1, 2 ... 6 */ /* has the same meaning as explained in the header comment li nes.*/ /* Parameter adjustments */ pdf -= -6; /* Function Body */ *irt = 0; for (i__ = -6; i__ <= 2; ++i__) { if (i__ <= 0) { pdf[i__] = ctq3df_(iset, &i__, x, q, &irt1); pdf[-i__] = pdf[i__]; } else { pdf[i__] = ctq3df_(iset, &i__, x, q, &irt1) + pdf[-i__]; } *irt += irt1; /* L10: */ } return 0; /* ************************* */ } /* ctq3pds_ */ doublereal ctq3df_0_(n__, iset, iprtn, xx, qq, irt, iorder, neff) int n__; integer *iset, *iprtn; doublereal *xx, *qq; integer *irt, *iorder, *neff; { /* Initialized data */ static doublereal alm[3] = { .239,.177,.247 }; static doublereal vlm[9] /* was [3][3] */ = { .239,.158,.063,.177,.132, .066,.247,.164,.066 }; static doublereal qms[9] /* was [3][3] */ = { 1.6,5.,180.,1.6,5.,180., 1.6,5.,180. }; static doublereal qmn[3] = { 1.6,1.6,1.6 }; static integer ist = 0; static integer lp = -10; static doublereal qsto = 1.2345; static integer iord[3] = { 2,1,2 }; /* System generated locals */ doublereal ret_val, d__1, d__2; /* Builtin functions */ double log(), exp(), pow_dd(); /* Subroutine */ int s_stop(); /* Local variables */ static doublereal alam; static integer iflv; static doublereal x, a0, a1, a2, a3, a4, a5, sb; static integer ip; static doublereal qi, sb2, sb3, sbl; /* Returns xf(x,Q) -- the momentum fraction distribution !! */ /* Returns valence and sea rather than combined flavor distr. */ /* Iset : PDF set label */ /* Iprtn : Parton label: 2, 1 = d_ and u_ valence */ /* 0 = gluon */ /* -1, ... -6 = u, d, s, c, b, t sea quarks */ /* XX : Bjorken-x */ /* QQ : scale parameter "Q" */ /* Irt : Return code */ /* 0 : no error */ /* 1 : parametrization is slightly negative; reset to 0.0. */ /* (This condition happens rarely -- only for large x where the */ /* absolute value of the parton distribution is extremely small.) */ /* --------- CTEQ3M */ switch(n__) { case 1: goto L_wlamd3; } /* --------- CTEQ3L */ /* --------- CTEQ3D */ x = *xx; *irt = 0; if (*iset == ist && qsto == *qq) { /* if only change is in x: */ if (*iprtn == lp) { goto L100; } /* if change in flv is within "light" partons: */ if (*iprtn >= -3 && lp >= -3) { goto L501; } } ip = abs(*iprtn); /* Set up Qi for SB */ if (ip >= 4) { if (*qq <= qms[ip + *iset * 3 - 7]) { ret_val = (float)0.; return ret_val; } qi = qms[ip + *iset * 3 - 7]; } else { qi = qmn[*iset - 1]; } /* Use "standard lambda" of parametrization program */ alam = alm[*iset - 1]; sbl = log(*qq / alam) / log(qi / alam); sb = log(sbl); sb2 = sb * sb; sb3 = sb2 * sb; L501: iflv = 3 - *iprtn; switch ((int)*iset) { case 1: goto L1; case 2: goto L2; case 3: goto L3; case 4: goto L311; } L1: switch ((int)iflv) { case 1: goto L11; case 2: goto L12; case 3: goto L13; case 4: goto L14; case 5: goto L15; case 6: goto L16; case 7: goto L17; case 8: goto L18; case 9: goto L19; } /* Ifl = 2 */ L11: a0 = exp((float)-.7266 - sb * (float)1.584 + sb2 * (float)1.259 - sb3 * ( float).04305); a1 = (float).5285 - sb * (float).3721 + sb2 * (float).515 - sb3 * (float) .1697; a2 = sb * (float).8282 + (float)4.075 - sb2 * (float).4496 + sb3 * (float) .2107; a3 = sb * (float)5.066 + (float)3.279 - sb2 * (float)9.134 + sb3 * (float) 2.897; a4 = (float).4399 - sb * (float).5888 + sb2 * (float).4802 - sb3 * (float) .1664; a5 = (float).3678 - sb * (float).8929 + sb2 * (float)1.592 - sb3 * (float) .5713; goto L100; /* Ifl = 1 */ L12: a0 = exp(sb * (float).1237 + (float).2259 + sb2 * (float).3035 - sb3 * ( float).2935); a1 = sb * (float).01651 + (float).5085 - sb2 * (float).03592 + sb3 * ( float).02782; a2 = sb * (float).4901 + (float)3.732 + sb2 * (float).2218 - sb3 * (float) .1116; a3 = (float)7.011 - sb * (float)6.62 + sb2 * (float)2.557 - sb3 * (float) .136; a4 = (float).8969 - sb * (float).2429 + sb2 * (float).1811 - sb3 * (float) .06888; a5 = sb * (float).2558 + (float).08636 - sb2 * (float).3082 + sb3 * ( float).2535; goto L100; /* Ifl = 0 */ L13: a0 = exp((float)-.2318 - sb * (float).9779 - sb2 * (float).3783 + sb3 * ( float).01037); a1 = sb * (float).1754 - (float).2916 - sb2 * (float).1884 + sb3 * (float) .06116; a2 = sb * (float).746 + (float)5.349 + sb2 * (float).2319 - sb3 * (float) .2622; a3 = (float)6.92 - sb * (float)3.454 + sb2 * (float)2.027 - sb3 * (float) .7626; a4 = sb * (float).1423 + (float)1.013 - sb2 * (float).1798 + sb3 * (float) .01872; a5 = sb * (float)2.303 - (float).05465 - sb2 * (float).9584 + sb3 * ( float).3098; goto L100; /* Ifl = -1 */ L14: a0 = exp((float)-2.906 - sb * (float).1069 - sb2 * (float)1.055 + sb3 * ( float).2496); a1 = sb * (float).06571 - (float).2875 - sb2 * (float).01987 - sb3 * ( float).0018; a2 = (float)9.854 - sb * (float).2715 - sb2 * (float).7407 + sb3 * (float) .2888; a3 = (float)15.83 - sb * (float)7.687 + sb2 * (float)3.428 - sb3 * (float) .3327; a4 = sb * (float).07599 + (float).9763 - sb2 * (float).2128 + sb3 * ( float).06852; a5 = sb * (float).9434 - (float).008444 + sb2 * (float).4152 - sb3 * ( float).1481; goto L100; /* Ifl = -2 */ L15: a0 = exp((float)-2.328 - sb * (float)3.061 + sb2 * (float)3.62 - sb3 * ( float)1.602); a1 = sb * (float).3198 - (float).3358 - sb2 * (float).421 + sb3 * (float) .1571; a2 = (float)8.478 - sb * (float)3.112 + sb2 * (float)5.243 - sb3 * (float) 2.255; a3 = sb * (float).3389 + (float)19.71 - sb2 * (float)5.268 + sb3 * (float) 2.099; a4 = (float)1.128 - sb * (float).4701 + sb2 * (float).7779 - sb3 * (float) .3506; a5 = sb * (float)3.341 - (float).4708 - sb2 * (float)3.375 + sb3 * (float) 1.353; goto L100; /* Ifl = -3 */ L16: a0 = exp(sb * (float)2.499 - (float)3.78 - sb2 * (float)4.962 + sb3 * ( float)1.936); a1 = (float)-.2639 - sb * (float).1575 + sb2 * (float).3584 - sb3 * ( float).1646; a2 = sb * (float)2.794 + (float)8.082 - sb2 * (float)5.438 + sb3 * (float) 2.321; a3 = (float)18.11 - sb * (float)20. + sb2 * (float)19.51 - sb3 * (float) 6.904; a4 = sb * (float).4972 + (float).9822 - sb2 * (float).869 + sb3 * (float) .3415; a5 = (float).1772 - sb * (float).6078 + sb2 * (float)3.341 - sb3 * (float) 1.473; goto L100; /* Ifl = -4 */ L17: a0 = pow_dd(&sb, &c_b20) * exp((float)-4.232 - sb * (float)1.808 + sb2 * ( float).5348); a1 = sb * (float).5846 - (float).2824 - sb2 * (float).723 + sb3 * (float) .2419; a2 = (float)5.683 - sb * (float)2.948 + sb2 * (float)5.916 - sb3 * (float) 2.56; a3 = sb * (float)4.795 + (float)2.051 - sb2 * (float)4.271 + sb3 * (float) .4174; a4 = sb * (float)1.717 + (float).1737 - sb2 * (float)1.978 + sb3 * (float) .6643; a5 = sb * (float)3.5 + (float).8689 - sb2 * (float)3.283 + sb3 * (float) 1.026; goto L100; /* Ifl = -5 */ L18: a0 = pow_dd(&sb, &c_b21) * exp((float)-1.496 - sb * (float)6.576 + sb2 * ( float)1.569); a1 = (float)-.214 - sb * (float).06419 - sb2 * (float).002741 + sb3 * ( float).003185; a2 = sb * (float).1049 + (float)5.781 - sb2 * (float).393 + sb3 * (float) .5174; a3 = sb * (float).5511 - (float).942 + sb2 * (float).8817 + sb3 * (float) 1.903; a4 = sb * (float).04232 + (float).02418 - sb2 * (float).01244 - sb3 * ( float).02365; a5 = sb * (float)1.794 + (float).7664 - sb2 * (float).4917 - sb3 * (float) .1284; goto L100; /* Ifl = -6 */ L19: a0 = pow_dd(&sb, &c_b22) * exp(sb * (float)1.154 - (float)8.46 + sb2 * ( float)8.838); a1 = (float)-.04316 - sb * (float).2976 + sb2 * (float).3174 - sb3 * ( float)1.429; a2 = sb * (float)2.273 + (float)4.91 + sb2 * (float)5.631 - sb3 * (float) 19.94; a3 = (float)11.9 - sb * (float)20. - sb2 * (float)20. + sb3 * (float) 12.92; a4 = (float).5771 - sb * (float).2552 + sb2 * (float).751 + sb3 * (float) .6923; a5 = (float)4.402 - sb * (float)1.627 - sb2 * (float)2.085 - sb3 * (float) 6.737; goto L100; L2: switch ((int)iflv) { case 1: goto L21; case 2: goto L22; case 3: goto L23; case 4: goto L24; case 5: goto L25; case 6: goto L26; case 7: goto L27; case 8: goto L28; case 9: goto L29; } /* Ifl = 2 */ L21: a0 = exp(sb * (float).4764 + (float).1141 - sb2 * (float)1.745 + sb3 * ( float).7728); a1 = (float).4275 - sb * (float).129 + sb2 * (float).3609 - sb3 * (float) .1689; a2 = sb * (float)2.946 + (float)3. - sb2 * (float)4.117 + sb3 * (float) 1.989; a3 = sb * (float)2.322 - (float)1.302 - sb2 * (float)4.258 + sb3 * (float) 2.109; a4 = (float)2.586 - sb * (float).192 - sb2 * (float).3754 + sb3 * (float) .2731; a5 = (float)-.2251 - sb * (float).5374 + sb2 * (float)2.245 - sb3 * ( float)1.034; goto L100; /* Ifl = 1 */ L22: a0 = exp(sb * (float).04205 + (float).1907 + sb2 * (float).2752 - sb3 * ( float).3171); a1 = sb * (float).02331 + (float).4611 - sb2 * (float).03403 + sb3 * ( float).03174; a2 = sb * (float).5739 + (float)3.504 + sb2 * (float).2676 - sb3 * (float) .1553; a3 = (float)7.452 - sb * (float)6.742 + sb2 * (float)2.849 - sb3 * (float) .1964; a4 = (float)1.116 - sb * (float).3435 + sb2 * (float).2865 - sb3 * (float) .1288; a5 = sb * (float).2714 + (float).06659 - sb2 * (float).2688 + sb3 * ( float).2763; goto L100; /* Ifl = 0 */ L23: a0 = exp((float)-.7631 - sb * (float).7241 - sb2 * (float)1.17 + sb3 * ( float).5343); a1 = sb * (float).3469 - (float).3573 - sb2 * (float).3396 + sb3 * (float) .09188; a2 = sb * (float).7458 + (float)5.604 - sb2 * (float).5082 + sb3 * (float) .1844; a3 = (float)15.49 - sb * (float)18.09 + sb2 * (float)11.62 - sb3 * (float) 3.483; a4 = sb * (float).1364 + (float).9881 - sb2 * (float).4421 + sb3 * (float) .2051; a5 = sb * (float)3.259 - (float).09505 - sb2 * (float)1.547 + sb3 * ( float).2918; goto L100; /* Ifl = -1 */ L24: a0 = exp((float)-2.74 - sb * (float).07987 - sb2 * (float).9015 - sb3 * ( float).09872); a1 = sb * (float).1244 - (float).3909 - sb2 * (float).04487 + sb3 * ( float).01277; a2 = sb * (float).2823 + (float)9.163 - sb2 * (float).772 - sb3 * (float) .00936; a3 = (float)10.8 - sb * (float)3.915 - sb2 * (float)1.153 + sb3 * (float) 2.649; a4 = (float).9894 - sb * (float).1647 - sb2 * (float).009426 + sb3 * ( float).002945; a5 = sb * (float).6998 - (float).3395 + sb2 * (float).7 - sb3 * (float) .0673; goto L100; /* Ifl = -2 */ L25: a0 = exp((float)-2.449 - sb * (float)3.513 + sb2 * (float)4.529 - sb3 * ( float)2.031); a1 = sb * (float).3411 - (float).405 - sb2 * (float).3669 + sb3 * (float) .1109; a2 = (float)7.47 - sb * (float)2.982 + sb2 * (float)5.503 - sb3 * (float) 2.419; a3 = sb * (float)1.638 + (float)15.03 - sb2 * (float)8.772 + sb3 * (float) 3.852; a4 = (float)1.137 - sb * (float)1.006 + sb2 * (float)1.485 - sb3 * (float) .6389; a5 = sb * (float)3.16 - (float).5299 - sb2 * (float)3.104 + sb3 * (float) 1.219; goto L100; /* Ifl = -3 */ L26: a0 = exp(sb * (float)1.25 - (float)3.64 - sb2 * (float)2.914 + sb3 * ( float).839); a1 = (float)-.3595 - sb * (float).05259 + sb2 * (float).3122 - sb3 * ( float).1642; a2 = sb * (float).9727 + (float)7.305 - sb2 * (float).9788 - sb3 * (float) .05193; a3 = (float)11.98 - sb * (float)17.99 + sb2 * (float)26.14 - sb3 * (float) 10.91; a4 = (float).9882 - sb * (float).6101 + sb2 * (float).9737 - sb3 * (float) .4935; a5 = (float)-.1186 - sb * (float).3231 + sb2 * (float)3.074 - sb3 * ( float)1.274; goto L100; /* Ifl = -4 */ L27: a0 = pow_dd(&sb, &c_b20) * exp((float)-3.718 - sb * (float)1.335 + sb2 * ( float).01651); a1 = sb * (float).7509 - (float).4719 - sb2 * (float).842 + sb3 * (float) .2901; a2 = (float)6.194 - sb * (float)1.641 + sb2 * (float)4.907 - sb3 * (float) 2.523; a3 = (float)4.426 - sb * (float)4.27 + sb2 * (float)6.581 - sb3 * (float) 3.474; a4 = sb * (float).9876 + (float).2683 - sb2 * (float).7612 + sb3 * (float) .178; a5 = sb * (float)4.41 - (float).4547 - sb2 * (float)3.712 + sb3 * (float) 1.245; goto L100; /* Ifl = -5 */ L28: a0 = pow_dd(&sb, &c_b33) * exp((float)-2.548 - sb * (float)7.66 + sb2 * ( float)3.702); a1 = (float)-.3122 - sb * (float).212 + sb2 * (float).5716 - sb3 * (float) .3773; a2 = (float)6.257 - sb * (float).08214 - sb2 * (float)2.537 + sb3 * ( float)2.981; a3 = sb * (float)2.131 - (float).6723 + sb2 * (float)9.599 - sb3 * (float) 7.91; a4 = sb * (float).04295 + (float).09169 - sb2 * (float).5017 + sb3 * ( float).3811; a5 = sb * (float)2.656 + (float).2402 - sb2 * (float)1.586 + sb3 * (float) .288; goto L100; /* Ifl = -6 */ L29: a0 = pow_dd(&sb, &c_b34) * exp(sb * (float)3.05 - (float)6.934 - sb2 * ( float).6943); a1 = (float)-.1713 - sb * (float).5167 + sb2 * (float)1.241 - sb3 * ( float)1.703; a2 = sb * (float)3.023 + (float)6.169 - sb2 * (float)19.72 + sb3 * (float) 10.69; a3 = (float)4.439 - sb * (float)17.46 + sb2 * (float)12.25 + sb3 * (float) .835; a4 = (float).5458 - sb * (float).4586 + sb2 * (float).9089 - sb3 * (float) .4049; a5 = (float)3.207 - sb * (float)3.362 + sb2 * (float)5.877 - sb3 * (float) 7.659; goto L100; L3: switch ((int)iflv) { case 1: goto L31; case 2: goto L32; case 3: goto L33; case 4: goto L34; case 5: goto L35; case 6: goto L36; case 7: goto L37; case 8: goto L38; case 9: goto L39; } /* Ifl = 2 */ L31: a0 = exp(sb * (float).4914 + (float).3961 - sb2 * (float)1.728 + sb3 * ( float).7257); a1 = (float).4162 - sb * (float).1419 + sb2 * (float).368 - sb3 * (float) .1618; a2 = sb * (float)3.028 + (float)3.248 - sb2 * (float)4.307 + sb3 * (float) 1.92; a3 = sb * (float)2.184 - (float)1.1 - sb2 * (float)3.82 + sb3 * (float) 1.717; a4 = (float)2.082 - sb * (float).2756 + sb2 * (float).3043 - sb3 * (float) .126; a5 = (float)-.4822 - sb * (float).5706 + sb2 * (float)2.243 - sb3 * ( float).976; goto L100; /* Ifl = 1 */ L32: a0 = exp(sb * (float).05814 + (float).2148 + sb2 * (float).2734 - sb3 * ( float).2902); a1 = sb * (float).01657 + (float).481 - sb2 * (float).038 + sb3 * (float) .03125; a2 = sb * (float).3923 + (float)3.509 + sb2 * (float).401 - sb3 * (float) .1932; a3 = (float)7.055 - sb * (float)6.552 + sb2 * (float)3.466 - sb3 * (float) .5657; a4 = (float)1.061 - sb * (float).3453 + sb2 * (float).4089 - sb3 * (float) .1817; a5 = sb * (float).2548 + (float).08687 - sb2 * (float).2967 + sb3 * ( float).2647; goto L100; /* Ifl = 0 */ L33: a0 = exp((float)-.4665 - sb * (float).7554 - sb2 * (float).3323 - sb3 * ( float)2.734e-5); a1 = sb * (float).2395 - (float).3359 - sb2 * (float).2377 + sb3 * (float) .07059; a2 = sb * (float).6086 + (float)5.451 + sb2 * (float).08606 - sb3 * ( float).1425; a3 = (float)10.26 - sb * (float)9.352 + sb2 * (float)4.879 - sb3 * (float) 1.15; a4 = (float).9935 - sb * (float).05017 - sb2 * (float).01707 - sb3 * ( float).001464; a5 = sb * (float)2.305 - (float).0416 - sb2 * (float)1.063 + sb3 * (float) .3211; goto L100; /* Ifl = -1 */ L34: a0 = exp(sb * (float).2296 - (float)3.323 - sb2 * (float)1.109 + sb3 * ( float).2223); a1 = sb * (float).08847 - (float).341 - sb2 * (float).01111 - sb3 * ( float).005927; a2 = (float)9.753 - sb * (float).5182 - sb2 * (float).467 + sb3 * (float) .1921; a3 = (float)19.77 - sb * (float)16. + sb2 * (float)9.481 - sb3 * (float) 1.864; a4 = sb * (float).002839 + (float).9818 - sb2 * (float).1188 + sb3 * ( float).03584; a5 = sb * (float)1.004 - (float).07934 + sb2 * (float).3704 - sb3 * ( float).122; goto L100; /* Ifl = -2 */ L35: a0 = exp((float)-2.714 - sb * (float)2.868 + sb2 * (float)3.7 - sb3 * ( float)1.671); a1 = sb * (float).3341 - (float).3893 - sb2 * (float).3897 + sb3 * (float) .142; a2 = (float)8.359 - sb * (float)3.267 + sb2 * (float)5.327 - sb3 * (float) 2.245; a3 = (float)23.59 - sb * (float)5.669 - sb2 * (float)4.602 + sb3 * (float) 3.153; a4 = (float)1.106 - sb * (float).4745 + sb2 * (float).7739 - sb3 * (float) .3417; a5 = sb * (float)3.433 - (float).5557 - sb2 * (float)3.39 + sb3 * (float) 1.354; goto L100; /* Ifl = -3 */ L36: a0 = exp(sb * (float)2.855 - (float)3.985 - sb2 * (float)5.208 + sb3 * ( float)1.937); a1 = (float)-.3337 - sb * (float).115 + sb2 * (float).3691 - sb3 * (float) .1709; a2 = sb * (float)3.641 + (float)7.968 - sb2 * (float)6.599 + sb3 * (float) 2.642; a3 = (float)18.73 - sb * (float)19.99 + sb2 * (float)17.34 - sb3 * (float) 5.813; a4 = sb * (float).5082 + (float).9731 - sb2 * (float).878 + sb3 * (float) .3231; a5 = (float)-.05542 - sb * (float).4189 + sb2 * (float)3.309 - sb3 * ( float)1.439; goto L100; /* Ifl = -4 */ L37: a0 = pow_dd(&sb, &c_b44) * exp((float)-3.952 - sb * (float)1.901 + sb2 * ( float).5137); a1 = sb * (float).6055 - (float).3543 - sb2 * (float).6941 + sb3 * (float) .2278; a2 = (float)5.955 - sb * (float)2.629 + sb2 * (float)5.337 - sb3 * (float) 2.3; a3 = sb * (float)4.882 + (float)1.933 - sb2 * (float)3.81 + sb3 * (float) .229; a4 = sb * (float)1.655 + (float).1806 - sb2 * (float)1.893 + sb3 * (float) .6395; a5 = sb * (float)3.612 + (float).479 - sb2 * (float)3.152 + sb3 * (float) .9684; goto L100; /* Ifl = -5 */ L38: a0 = pow_dd(&sb, &c_b45) * exp((float)-1.825 - sb * (float)7.464 + sb2 * ( float)2.143); a1 = (float)-.2604 - sb * (float).14 + sb2 * (float).1702 - sb3 * (float) .08476; a2 = sb * (float).6275 + (float)6.005 - sb2 * (float)2.535 + sb3 * (float) 2.219; a3 = sb * (float)1.149 - (float).9067 + sb2 * (float)1.974 + sb3 * (float) 4.716; a4 = sb * (float).05945 + (float).03915 - sb2 * (float).09844 + sb3 * ( float).02783; a5 = sb * (float)1.994 + (float).55 - sb2 * (float).6727 - sb3 * (float) .151; goto L100; /* Ifl = -6 */ L39: a0 = pow_dd(&sb, &c_b46) * exp(sb * (float).3793 - (float)8.553 + sb2 * ( float)9.998); a1 = (float)-.0587 - sb * (float).2792 + sb2 * (float).6526 - sb3 * ( float)1.984; a2 = sb * (float).4473 + (float)4.716 + sb2 * (float)11.28 - sb3 * (float) 19.37; a3 = (float)12.89 - sb * (float)17.42 - sb2 * (float)19.83 - sb3 * (float) .9274; a4 = (float).5647 - sb * (float).2732 + sb2 * (float)1.074 + sb3 * (float) .5981; a5 = (float)4.39 - sb * (float)1.262 - sb2 * (float).9026 - sb3 * (float) 9.394; goto L100; L311: s_stop("This option is not currently supported.", 39L); L100: d__1 = 1. - x; d__2 = log(1. / x + 1.); ret_val = a0 * pow_dd(&x, &a1) * pow_dd(&d__1, &a2) * (a3 * pow_dd(&x, & a4) + 1.) * pow_dd(&d__2, &a5); if (ret_val < 0.) { ret_val = 0.; *irt = 1; } ist = *iset; lp = *iprtn; qsto = *qq; return ret_val; /* ----------------------- */ L_wlamd3: /* Returns the EFFECTIVE QCD lambda values for order=Iorder and */ /* effective # of flavors = Neff for each of the PDF sets. */ *iorder = iord[*iset - 1]; ret_val = vlm[*neff + *iset * 3 - 7]; return ret_val; /* ************************* */ } /* ctq3df_ */ doublereal ctq3df_(iset, iprtn, xx, qq, irt) integer *iset, *iprtn; doublereal *xx, *qq; integer *irt; { return ctq3df_0_(0, iset, iprtn, xx, qq, irt, (integer *)0, (integer *)0); } doublereal wlamd3_(iset, iorder, neff) integer *iset, *iorder, *neff; { return ctq3df_0_(1, iset, (integer *)0, (doublereal *)0, (doublereal *)0, (integer *)0, iorder, neff); }