Changeset 19917 for trunk/Mars

Timestamp:

12/17/19 17:26:49 (5 years ago)

Author:

tbretz

Message:

Added a possibility to produce a pre-defined wavelength distirbution.

Location:

trunk/Mars

Files:

• mbase/MSpline3.cc (modified) (1 diff)
• mbase/MSpline3.h (modified) (1 diff)
• msimreflector/MSimRays.cc (modified) (4 diffs)
• msimreflector/MSimRays.h (modified) (4 diffs)

trunk/Mars/mbase/MSpline3.cc

@@ -234 +234 @@
     }
 
-    return TMath::TwoPi()*MSpline3(x.GetArray(), y.GetArray(), n).Integral();
-}
+    return MSpline3(x.GetArray(), y.GetArray(), n).Integral();
+}
+
+MSpline3 MSpline3::GetIntegralSpline() const
+{
+    const Int_t n = GetNp();
+
+    MArrayD x(n+2);
+    MArrayD y(n+2);
+
+    GetKnot(0, x[0], y[0]);
+
+    y[0] = 0;
+
+    GetKnot(0, x[1], y[1]);
+
+    y[1] = 0;
+
+    for (int i=0; i<n-1; i++)
+    {
+        GetKnot(i+1, x[i+2], y[i+2]);
+        y[i+2] = y[i+1]+IntegralBin(i, x[i+2]);
+    }
+
+    for (int i=2; i<=n; i++)
+        y[i] /= y[n];
+
+    y[n+1] = 1;
+
+    x[0]   = x[1]  -(x[2]-x[1]);
+    x[n+1] = x[n]  +(x[n]-x[n-1]);
+
+    return MSpline3(x.GetArray(), y.GetArray(), n+2, "b1 e1");
+}
+
 
 

trunk/Mars/mbase/MSpline3.h

@@ -64 +64 @@
     Double_t IntegralSolidAngle() const;
 
+    MSpline3 GetIntegralSpline() const;
+
     void Scale(double=1);
 

trunk/Mars/msimreflector/MSimRays.cc

@@ -52 +52 @@
 #include "MParList.h"
 
+#include "MSpline3.h"
 #include "MQuaternion.h"
 
@@ -71 +72 @@
     : fEvt(0), fReflector(0), fPointPos(0), fSource(0),
     fNumPhotons(1000), fHeight(-1), fWavelengthMin(-1), fWavelengthMax(-1),
-    fNameReflector("MReflector"), fNamePointPos("MPointingPos"),
+    fRandomDist(0), fNameReflector("MReflector"), fNamePointPos("MPointingPos"),
     fNameSource("Source")
 {
     fName  = name  ? name  : "MSimRays";
     fTitle = title ? title : "Task to calculate reflection os a mirror";
+}
+
+MSimRays::~MSimRays()
+{
+    delete fRandomDist;
 }
 
@@ -234 +240 @@
         dat.SetDirection(w);
 
-        if (fWavelengthMin>0 && fWavelengthMax>0)
-            dat.SimWavelength(fWavelengthMin, fWavelengthMax);
+        if (fRandomDist)
+        {
+            dat.SetWavelength(fRandomDist->FindX(gRandom->Uniform()));
+        }
+        else
+        {
+            if (fWavelengthMin>0 && fWavelengthMax>0)
+                dat.SimWavelength(fWavelengthMin, fWavelengthMax);
+        }
 
         idx++;
@@ -247 +260 @@
 }
 
+void MSimRays::SetDistribution(const MSpline3 &s)
+{
+    delete fRandomDist;
+    fRandomDist = new MSpline3(s.GetIntegralSpline());
+}
+
+bool MSimRays::ReadDistribution(const char *filename, const char *fmt)
+{
+    TGraph g(filename, fmt);
+    if (g.GetN()<2)
+        return false;
+
+    SetDistribution(MSpline3(g));
+    return true;
+}
+
 // --------------------------------------------------------------------------
 //

trunk/Mars/msimreflector/MSimRays.h

@@ -7 +7 @@
 
 class MParList;
+class MSpline3;
 class MGeomCam;
 class MPointingPos;
@@ -28 +29 @@
     Double_t fWavelengthMax;
 
+    MSpline3 *fRandomDist;      //!
+
     TString fNameReflector;
     TString fNamePointPos;
@@ -41 +44 @@
 public:
     MSimRays(const char *name=NULL, const char *title=NULL);
+    ~MSimRays();
 
     void SetNameReflector(const char *name="MReflector")  { fNameReflector = name; }
@@ -50 +54 @@
 
     void SetWavelengthRange(Double_t min, Double_t max) { fWavelengthMin = min; fWavelengthMax = max; }
+    void SetDistribution(const MSpline3 &s);
+    bool ReadDistribution(const char *filename, const char *fmt="%lg %lg");
 
     // MSimRays