Changeset 1352 for trunk/WuerzburgSoft/Thomas/mphys/MElectron.cc

Timestamp:

06/10/02 08:49:28 (22 years ago)

Author:

tbretz

Message:

*** empty log message ***

File:

• trunk/WuerzburgSoft/Thomas/mphys/MElectron.cc (modified) (4 diffs)

trunk/WuerzburgSoft/Thomas/mphys/MElectron.cc

@@ -34 +34 @@
 #include <TF1.h>
 
+#include "MPhoton.h"
+
 ClassImp(MElectron);
 
@@ -238 +240 @@
 
 
-Double_t MElectron::EnergyLoss(Double_t *x, Double_t *k = NULL)
+Double_t MElectron::EnergyLoss(Double_t *x, Double_t *k, Double_t *ep)
 {
     Double_t E = x[0];
@@ -255 +257 @@
     Double_t e = pow(10, fP.GetRandom());
 
+    if (ep)
+        *ep = e;
+
     Double_t omega = e*E/E0/E0;
     Double_t Gamma = 4.*omega;
@@ -269 +274 @@
 }
 
-Double_t MElectron::GetEnergyLoss(Double_t E, Double_t z)
+Double_t MElectron::GetEnergyLoss(Double_t E, Double_t z, Double_t *ep)
 {
     return EnergyLoss(&E, &z);
 }
 
-Double_t MElectron::GetEnergyLoss() const
-{
-    return EnergyLoss((Double_t*)&fEnergy, (Double_t*)&fZ);
-}
+Double_t MElectron::GetEnergyLoss(Double_t *ep) const
+{
+    return EnergyLoss((Double_t*)&fEnergy, (Double_t*)&fZ, ep);
+}
+
+#include <iostream.h>
+
+MPhoton *MElectron::DoInvCompton()
+{
+    Double_t E0 = 511e-6; //[GeV]
+
+    Double_t epsilon;
+    Double_t e = GetEnergyLoss(&epsilon);
+    /*
+     Double_t gamma = fEnergy/E0;
+     Double_t beta = sqrt(1-1/(gamma*gamma));
+     Double_t bega = sqrt(gamma*gamma-1);
+
+     Double_t theta = TMath::Pi()/4;
+
+     // er: photon energy before interaction rest frame
+     // e:  photon energy after interaction lab
+     Double_t er = gamma*epsilon*(1-beta*cos(theta)); // photon energy rest frame
+     Double_t r1 = fEnergy/e;
+     Double_t r2 = E0/er;
+     Double_t ctheta = (r1-r2-1)/(beta*r1-1);
+     Double_t tg = sqrt(1-ctheta*ctheta)/gamma/(beta+ctheta);
+    */
+    fEnergy -= e;
+
+    MPhoton &p = *new MPhoton(e, fZ);
+    p = *this;
+
+    /*
+     cout << "(" << atan(tg)*180/TMath::Pi() << ")" << flush;
+
+     static TRandom rand(0);
+     p->SetNewDirection(atan(tg), rand.Uniform(TMath::Pi()*2));
+     */
+
+    // MISSING: Electron angle
+    //
+    // E1 = fEnergy  (after!)
+    //
+    // sin(t) = (epsilon sin(theta) - e sin(atan(tg)))/sqrt(E1*E1-E0*E0)
+
+    return &p;
+}