Changeset 19917


Ignore:
Timestamp:
Dec 17, 2019, 5:26:49 PM (3 months ago)
Author:
tbretz
Message:
Added a possibility to produce a pre-defined wavelength distirbution.
Location:
trunk/Mars
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/Mars/mbase/MSpline3.cc

    r19128 r19917  
    234234    }
    235235
    236     return TMath::TwoPi()*MSpline3(x.GetArray(), y.GetArray(), n).Integral();
    237 }
     236    return MSpline3(x.GetArray(), y.GetArray(), n).Integral();
     237}
     238
     239MSpline3 MSpline3::GetIntegralSpline() const
     240{
     241    const Int_t n = GetNp();
     242
     243    MArrayD x(n+2);
     244    MArrayD y(n+2);
     245
     246    GetKnot(0, x[0], y[0]);
     247
     248    y[0] = 0;
     249
     250    GetKnot(0, x[1], y[1]);
     251
     252    y[1] = 0;
     253
     254    for (int i=0; i<n-1; i++)
     255    {
     256        GetKnot(i+1, x[i+2], y[i+2]);
     257        y[i+2] = y[i+1]+IntegralBin(i, x[i+2]);
     258    }
     259
     260    for (int i=2; i<=n; i++)
     261        y[i] /= y[n];
     262
     263    y[n+1] = 1;
     264
     265    x[0]   = x[1]  -(x[2]-x[1]);
     266    x[n+1] = x[n]  +(x[n]-x[n-1]);
     267
     268    return MSpline3(x.GetArray(), y.GetArray(), n+2, "b1 e1");
     269}
     270
    238271
    239272
  • trunk/Mars/mbase/MSpline3.h

    r19128 r19917  
    6464    Double_t IntegralSolidAngle() const;
    6565
     66    MSpline3 GetIntegralSpline() const;
     67
    6668    void Scale(double=1);
    6769
  • trunk/Mars/msimreflector/MSimRays.cc

    r19788 r19917  
    5252#include "MParList.h"
    5353
     54#include "MSpline3.h"
    5455#include "MQuaternion.h"
    5556
     
    7172    : fEvt(0), fReflector(0), fPointPos(0), fSource(0),
    7273    fNumPhotons(1000), fHeight(-1), fWavelengthMin(-1), fWavelengthMax(-1),
    73     fNameReflector("MReflector"), fNamePointPos("MPointingPos"),
     74    fRandomDist(0), fNameReflector("MReflector"), fNamePointPos("MPointingPos"),
    7475    fNameSource("Source")
    7576{
    7677    fName  = name  ? name  : "MSimRays";
    7778    fTitle = title ? title : "Task to calculate reflection os a mirror";
     79}
     80
     81MSimRays::~MSimRays()
     82{
     83    delete fRandomDist;
    7884}
    7985
     
    234240        dat.SetDirection(w);
    235241
    236         if (fWavelengthMin>0 && fWavelengthMax>0)
    237             dat.SimWavelength(fWavelengthMin, fWavelengthMax);
     242        if (fRandomDist)
     243        {
     244            dat.SetWavelength(fRandomDist->FindX(gRandom->Uniform()));
     245        }
     246        else
     247        {
     248            if (fWavelengthMin>0 && fWavelengthMax>0)
     249                dat.SimWavelength(fWavelengthMin, fWavelengthMax);
     250        }
    238251
    239252        idx++;
     
    247260}
    248261
     262void MSimRays::SetDistribution(const MSpline3 &s)
     263{
     264    delete fRandomDist;
     265    fRandomDist = new MSpline3(s.GetIntegralSpline());
     266}
     267
     268bool MSimRays::ReadDistribution(const char *filename, const char *fmt)
     269{
     270    TGraph g(filename, fmt);
     271    if (g.GetN()<2)
     272        return false;
     273
     274    SetDistribution(MSpline3(g));
     275    return true;
     276}
     277
    249278// --------------------------------------------------------------------------
    250279//
  • trunk/Mars/msimreflector/MSimRays.h

    r19599 r19917  
    77
    88class MParList;
     9class MSpline3;
    910class MGeomCam;
    1011class MPointingPos;
     
    2829    Double_t fWavelengthMax;
    2930
     31    MSpline3 *fRandomDist;      //!
     32
    3033    TString fNameReflector;
    3134    TString fNamePointPos;
     
    4144public:
    4245    MSimRays(const char *name=NULL, const char *title=NULL);
     46    ~MSimRays();
    4347
    4448    void SetNameReflector(const char *name="MReflector")  { fNameReflector = name; }
     
    5054
    5155    void SetWavelengthRange(Double_t min, Double_t max) { fWavelengthMin = min; fWavelengthMax = max; }
     56    void SetDistribution(const MSpline3 &s);
     57    bool ReadDistribution(const char *filename, const char *fmt="%lg %lg");
    5258
    5359    // MSimRays
Note: See TracChangeset for help on using the changeset viewer.