106{
107
108
109
115 G4int chargeBalance =
G4lrint(projectile->GetDefinition()->GetPDGCharge()
116 + targets[0]->GetDefinition()->GetPDGCharge()
117 + targets[1]->GetDefinition()->GetPDGCharge());
118
119 G4int baryonBalance = projectile->GetDefinition()->GetBaryonNumber()
120 + targets[0]->GetDefinition()->GetBaryonNumber()
121 + targets[1]->GetDefinition()->GetBaryonNumber();
122
123
124
126 theT1 = toSPS * theT1;
127 theT2 = toSPS * theT2;
128 thePro = toSPS * thePro;
130
131
136 theT1 = toZ * theT1;
137 theT2 = toZ * theT2;
138 thePro = toZ * thePro;
140
141
145 {
147 temp=d1;d1=d2;d2=temp;
148 }
151 {
153 {
155 {
157 }
158 else
159 {
161 }
162 }
164 {
166 }
167 }
169 {
171 {
173 {
175 }
176 else
177 {
179 }
180 }
182 {
184 }
185 }
186
187
188 G4double M_sq = (thePro+theT1+theT2).mag2();
192 G4double p = std::sqrt((m_sq*m_sq - 4.*m1_sq * m2_sq)/(4.*M_sq));
195 G4ThreeVector pFinal(p*std::sin(std::acos(costh))*std::cos(phi), p*std::sin(std::acos(costh))*std::sin(phi), p*costh);
196
197
198
199 G4double eFinal1 = std::sqrt(m1_sq+pFinal.mag2());
201 G4double eFinal2 = std::sqrt(m2_sq+pFinal.mag2());
203
204
205 final1 = toLab * final1;
206 final2 = toLab * final2;
207
208
209 final1 = fromSPS * final1;
210 final2 = fromSPS * final2;
211
212
218 result->push_back(f1);
219 result->push_back(f2);
220
221 for(size_t hpw=0; hpw<result->size(); hpw++)
222 {
223 energyBalance-=result->operator[](hpw)->Get4Momentum().t();
224 chargeBalance-=
G4lrint(result->operator[](hpw)->GetDefinition()->GetPDGCharge());
225 baryonBalance-=result->operator[](hpw)->GetDefinition()->GetBaryonNumber();
226 }
227 if(std::getenv("AbsorptionEnergyBalanceCheck"))
228 std::cout << "DEBUGGING energy balance B: "
229 <<energyBalance<<" "
230 <<chargeBalance<<" "
231 <<baryonBalance<<" "
233 return result;
234}
HepLorentzRotation & rotateY(double delta)
HepLorentzRotation & rotateZ(double delta)
HepLorentzRotation inverse() const
static G4Neutron * NeutronDefinition()
G4double GetPDGMass() const
G4double GetPDGCharge() const
static G4Proton * ProtonDefinition()