56 {
57
59
60
64
67 <<_nA<<","<<_nB<<","<<_nC<<endl;
68 }
69
70
74
77 <<_JA2<<","<<_JB2<<","<<_JC2<<endl;
78 }
79
80
81
82 int* _lambdaA2=new int[_nA];
83 int* _lambdaB2=new int[_nB];
84 int* _lambdaC2=new int[_nC];
85
87 int ib,ic;
88 for(ib=0;ib<_nB;ib++){
90 }
91
92
93 int i;
94
95
96 fillHelicity(_lambdaA2,_nA,_JA2);
97 fillHelicity(_lambdaB2,_nB,_JB2);
98 fillHelicity(_lambdaC2,_nC,_JC2);
99
101 report(
INFO,
"EvtGen")<<
"Helicity states of particle A:"<<endl;
102 for(i=0;i<_nA;i++){
104 }
105
106 report(
INFO,
"EvtGen")<<
"Helicity states of particle B:"<<endl;
107 for(i=0;i<_nB;i++){
109 }
110
111 report(
INFO,
"EvtGen")<<
"Helicity states of particle C:"<<endl;
112 for(i=0;i<_nC;i++){
114 }
115
116 report(
INFO,
"EvtGen")<<
"Will now figure out the valid (M_LS) states:"<<endl;
117
118 }
119
120 int Lmin=std::max(_JA2-_JB2-_JC2,std::max(_JB2-_JA2-_JC2,_JC2-_JA2-_JB2));
121 if (Lmin<0) Lmin=0;
122
123 int Lmax=_JA2+_JB2+_JC2;
124
125 int L;
126
127 int _nPartialWaveAmp=0;
128
129 int _nL[50];
130 int _nS[50];
131
132 for (L=Lmin;L<=Lmax;L+=2){
133 int Smin=
abs(L-_JA2);
134 if (Smin<
abs(_JB2-_JC2)) Smin=
abs(_JB2-_JC2);
135 int Smax=L+_JA2;
136 if (Smax>
abs(_JB2+_JC2)) Smax=
abs(_JB2+_JC2);
137 int S;
138 for (S=Smin;S<=Smax;S+=2){
139 _nL[_nPartialWaveAmp]=L;
140 _nS[_nPartialWaveAmp]=S;
141
142 _nPartialWaveAmp++;
144 report(
INFO,
"EvtGen")<<
"M["<<L<<
"]["<<S<<
"]"<<endl;
145 }
146 }
147 }
148
150
151 int argcounter=0;
152
154
155 for(i=0;i<_nPartialWaveAmp;i++){
157 argcounter+=2;
159 report(
INFO,
"EvtGen")<<
"M["<<_nL[i]<<
"]["<<_nS[i]<<
"]="<<_M[i]<<endl;
160 }
161 }
162
163
164
165
166
167 for(ib=0;ib<_nB;ib++){
168 for(ic=0;ic<_nC;ic++){
169 _HBC[ib][ic]=0.0;
170 if (
abs(_lambdaB2[ib]-_lambdaC2[ic])<=_JA2){
171 for(i=0;i<_nPartialWaveAmp;i++){
172 int L=_nL[i];
173 int S=_nS[i];
174 int lambda2=_lambdaB2[ib];
175 int lambda3=_lambdaC2[ic];
176 int s1=_JA2;
177 int s2=_JB2;
178 int s3=_JC2;
179 int m1=lambda2-lambda3;
182
184 report(
INFO,
"EvtGen") <<
"s2,lambda2:"<<s2<<
" "<<lambda2<<endl;
185 }
186
187 double fkwTmp = (L+1.0)/(s1+1.0);
188
190
192 *c1.coef(S,m1,s2,s3,lambda2,-lambda3)
193 *c2.coef(s1,m1,L,S,0,m1)*_M[i];
194
195 _HBC[ib][ic]+=tmp;
196 }
197 }
199 report(
INFO,
"EvtGen")<<
"_HBC["<<ib<<
"]["<<ic<<
"]="<<_HBC[ib][ic]<<endl;
200 }
201 }
202 }
203 }
204
208 _HBC);
209
210
211}
EvtComplex exp(const EvtComplex &c)
ostream & report(Severity severity, const char *facility)
void checkNDaug(int d1, int d2=-1)
void checkNArg(int a1, int a2=-1, int a3=-1, int a4=-1)
static EvtSpinType::spintype getSpinType(EvtId i)
static int getSpin2(spintype stype)
static int getSpinStates(spintype stype)