87 using namespace CLHEP;
92 {2190., 1920., 1700., 1600., 1440., 1232. };
106 115.603, 133.424, 153.991, 177.729, 205.126, 236.746, 273.24, 315.361, 363.973, 420.08, 484.836, 559.573, 645.832,
107 745.387, 860.289, 992.903, 1145.96, 1322.61, 1526.49, 1761.8, 2033.38, 2346.83, 2708.59, 3126.12, 3608.02, 4164.19,
108 4806.1, 5546.97, 6402.04, 7388.91, 8527.92, 9842.5, 11359.7, 13110.8, 15131.9, 17464.5, 20156.6, 23263.8, 26849.9,
109 30988.8, 35765.7, 41279, 47642.2, 54986.3, 63462.4, 73245.2, 84536, 97567.2, 112607, 129966 };
117 #ifdef G4MULTITHREADED
137 outFile <<
"G4NuMuNucleusCcModel is a neutrino-nucleus (charge current) scattering\n"
138 <<
"model which uses the standard model \n"
139 <<
"transfer parameterization. The model is fully relativistic\n";
151 G4int nSize(0), i(0), j(0),
k(0);
155 #ifdef G4MULTITHREADED
161 #ifdef G4MULTITHREADED
169 char* path = getenv(
"G4PARTICLEXSDATA");
170 std::ostringstream ost1, ost2, ost3, ost4;
171 ost1 << path <<
"/" <<
"neutrino" <<
"/" << pName <<
"/xarraycckr";
173 std::ifstream filein1( ost1.str().c_str() );
179 for( k = 0; k <
fNbin; ++
k )
181 for( i = 0; i <=
fNbin; ++i )
189 ost2 << path <<
"/" <<
"neutrino" <<
"/" << pName <<
"/xdistrcckr";
190 std::ifstream filein2( ost2.str().c_str() );
194 for( k = 0; k <
fNbin; ++
k )
196 for( i = 0; i <
fNbin; ++i )
204 ost3 << path <<
"/" <<
"neutrino" <<
"/" << pName <<
"/q2arraycckr";
205 std::ifstream filein3( ost3.str().c_str() );
209 for( k = 0; k <
fNbin; ++
k )
211 for( i = 0; i <=
fNbin; ++i )
213 for( j = 0; j <=
fNbin; ++j )
222 ost4 << path <<
"/" <<
"neutrino" <<
"/" << pName <<
"/q2distrcckr";
223 std::ifstream filein4( ost4.str().c_str() );
227 for( k = 0; k <
fNbin; ++
k )
229 for( i = 0; i <=
fNbin; ++i )
231 for( j = 0; j <
fNbin; ++j )
303 G4double cost(1.), sint(0.),
phi(0.), muMom(0.), massX2(0.), massX(0.), massR(0.), eCut(0.);
309 G4int pdgP(0), qB(0);
321 sint = std::sqrt( (1.0 - cost)*(1.0 + cost) );
356 if( pName ==
"nu_mu" ) pdgP = 211;
361 eCut = (
fMpi + mTarg)*(
fMpi + mTarg) - (massX + massR)*(massX + massR);
365 if ( lvX.
e() > eCut )
383 sint = std::sqrt( (1.0 - cost)*(1.0 + cost) );
424 if( pName ==
"nu_mu" ) qB = 2;
458 if( pName ==
"nu_mu" )
476 if( pName ==
"nu_mu" )
500 if( pName ==
"nu_mu" )
528 else if ( eX < 95000.*
GeV )
530 if (
fProton && pName ==
"nu_mu" ) qB = 2;
531 else if(
fProton && pName ==
"anti_nu_mu" ) qB = 0;
532 else if( !
fProton && pName ==
"nu_mu" ) qB = 1;
533 else if( !
fProton && pName ==
"anti_nu_mu" ) qB = -1;
541 if( pName ==
"nu_mu" ) pdgP = 211;
571 G4double e3(0.), pMu2(0.), pX2(0.), nMom(0.), rM(0.), hM(0.), tM = targetNucleus.
AtomicMass(A,Z);
603 pMu2 = fEmu*fEmu -
fMu*
fMu;
612 if( iTer >= iTerMax ) {
fBreak =
true;
return; }
627 sint = std::sqrt( (1.0 - cost)*(1.0 + cost) );
646 ei = tM - sqrt( (rM + Ex)*(rM + Ex) + nMom*nMom );
649 nm2 =
ei*
ei - nMom*nMom;
652 while( nm2 < 0. && iTer < iTerMax );
654 if( iTer >= iTerMax ) {
fBreak =
true;
return; }
707 pMu2 = fEmu*fEmu -
fMu*
fMu;
721 if( iTer >= iTerMax ) {
fBreak =
true;
return; }
735 sint = std::sqrt( (1.0 - cost)*(1.0 + cost) );
753 G4int i(0), nBin(50);
756 for( i = 0; i < nBin; ++i )
781 else xx = x1 + (e-
e1)*(x2-x1)/(e2-
e1);
795 for( i = 0; i < nBin; ++i )
823 else xx = x1 + (prob-p1)*(x2-x1)/(p2-p1);
842 qq1 =
GetQkr( 0, jX, prob);
844 else if ( iE >= nBin-1)
846 qq1 =
GetQkr( nBin-1, jX, prob);
858 else qq1 = q1 + (e-
e1)*(q2-q1)/(e2-
e1);
865 qq2 =
GetQkr( iE, 0, prob);
867 else if ( jX >= nBin)
869 qq2 =
GetQkr( iE, nBin, prob);
881 else qq2 = q1 + (e-
e1)*(q2-q1)/(e2-
e1);
897 for( i = 0; i < nBin; ++i )
925 else qq = q1 + (prob-p1)*(q2-q1)/(p2-p1);