Ignore:
Timestamp:
08/31/04 10:21:26 (20 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mimage
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/mimage/MHHillasSrc.cc

    r2438 r4817  
    7070    // connect all the histogram with the container fHist
    7171    //
    72     fAlpha = new TH1F("Alpha", "Alpha of Ellipse",            181, -90,  90);
    73     fDist  = new TH1F("Dist",  "Dist of Ellipse",             100,   0, 445);
    74     fCosDA = new TH1F("CosDA", "cos(Delta,Alpha) of Ellipse", 101,  -1,   1);
     72    fAlpha = new TH1F("Alpha", "Alpha of Ellipse",                181, -90,  90);
     73    fDist  = new TH1F("Dist",  "Dist of Ellipse",                 100,   0, 445);
     74    fCosDA = new TH1F("CosDA", "cos(Delta,Alpha) of Ellipse",     101,  -1,   1);
     75    fDCA   = new TH1F("DCA",   "Distance of closest aproach",     101,  -1,   1);
     76    fDCADelta  = new TH1F("DCADelta",  "Angle between shower and x-axis", 101,   0, 360);
    7577
    7678    fAlpha->SetDirectory(NULL);
    7779    fDist->SetDirectory(NULL);
    7880    fCosDA->SetDirectory(NULL);
     81    fDCA->SetDirectory(NULL);
     82    fDCADelta->SetDirectory(NULL);
    7983
    8084    fAlpha->SetXTitle("\\alpha [\\circ]");
    8185    fDist->SetXTitle("Dist [mm]");
    8286    fCosDA->SetXTitle("cos(\\delta,\\alpha)");
     87    fDCA->SetXTitle("DCA [\\circ]");
     88    fDCADelta->SetXTitle("DCADelta [0, 2\\pi]");
    8389
    8490    fAlpha->SetYTitle("Counts");
    8591    fDist->SetYTitle("Counts");
    8692    fCosDA->SetYTitle("Counts");
     93    fDCA->SetYTitle("Counts");
     94    fDCADelta->SetYTitle("Counts");
    8795}
    8896
     
    96104    delete fDist;
    97105    delete fCosDA;
     106    delete fDCA;
     107    delete fDCADelta;
    98108}
    99109
     
    121131    ApplyBinning(*plist, "Alpha",    fAlpha);
    122132    ApplyBinning(*plist, "Dist",     fDist);
     133    ApplyBinning(*plist, "DCA",      fDCA);
     134    ApplyBinning(*plist, "DCADelta",     fDCADelta);
    123135
    124136    return kTRUE;
     
    143155    fDist ->Fill(fUseMmScale ? h.GetDist() : fMm2Deg*h.GetDist(), w);
    144156    fCosDA->Fill(h.GetCosDeltaAlpha(), w);
     157    fDCA  ->Fill(fUseMmScale ? h.GetDCA() : fMm2Deg*h.GetDCA(), w);
     158    fDCADelta ->Fill(h.GetDCADelta(), w);
    145159
    146160    return kTRUE;
     
    186200    const Double_t scale = mmscale ? 1./fMm2Deg : fMm2Deg;
    187201    MH::ScaleAxis(fDist, scale);
     202    MH::ScaleAxis(fDCA,  scale);
    188203
    189204    fDist->SetXTitle(mmscale ? "Dist [mm]" : "Dist [\\circ]");
     205    fDCA->SetXTitle(mmscale ? "DCA [mm]" : "DCA [\\circ]");
    190206
    191207    fUseMmScale = mmscale;
     
    217233    fDist->Draw();
    218234
    219     delete pad->GetPad(3);
     235    pad->cd(3);
     236    gPad->SetBorderMode(0);
     237    fDCA->Draw();
    220238
    221239    pad->cd(4);
    222240    gPad->SetBorderMode(0);
    223     //gPad->SetLogy();
     241
     242    TVirtualPad *p = gPad;
     243    p->Divide(1,2);
     244    p->cd(1);
     245    gPad->SetBorderMode(0);
    224246    fCosDA->Draw();
     247
     248    p->cd(2);
     249    gPad->SetBorderMode(0);
     250    fDCADelta->Draw();
    225251
    226252    pad->Modified();
  • trunk/MagicSoft/Mars/mimage/MHHillasSrc.h

    r2416 r4817  
    1515    TH1F *fDist;      //->
    1616    TH1F *fCosDA;     //->
     17
     18    TH1F *fDCA;       //->
     19    TH1F *fDCADelta;  //->
    1720
    1821    Float_t fMm2Deg;
     
    3437    TH1F *GetHistDist()          { return fDist; }
    3538    TH1F *GetHistCosDeltaAlpha() { return fCosDA; }
     39    TH1F *GetHistDCA()           { return fDCA; }
     40    TH1F *GetHistDCADelta()      { return fDCADelta; }
    3641
    3742    void Draw(Option_t *opt=NULL);
  • trunk/MagicSoft/Mars/mimage/MHillasSrc.cc

    r4710 r4817  
    5757// Version 4:
    5858// ----------
    59 //
    6059// fHeadTail        removed
     60//
     61//
     62// Version 5:
     63// ----------
     64//  - added Float_t fDCA;      // [mm]   Distance to closest approach 'DCA'
     65//  - added Float_t fDCADelta; // [deg]  Angle of the shower axis with respect
     66//                                       to the x-axis [0,2pi]
    6167//
    6268/////////////////////////////////////////////////////////////////////////////
     
    9197    fAlpha         =  0;
    9298    fCosDeltaAlpha =  0;
     99
     100    fDCA           = -1;
     101    fDCADelta          =  0;
    93102}
    94103
     
    99108//  you call the Reset member function before.
    100109//
    101 Bool_t MHillasSrc::Calc(const MHillas &hillas)
    102 {
    103     const Double_t mx = hillas.GetMeanX();     // [mm]
    104     const Double_t my = hillas.GetMeanY();     // [mm]
    105 
    106     const Double_t sx = mx - fSrcPos->GetX();   // [mm]
    107     const Double_t sy = my - fSrcPos->GetY();   // [mm]
    108 
    109     const Double_t sd = hillas.GetSinDelta();  // [1]
    110     const Double_t cd = hillas.GetCosDelta();  // [1]
    111 
    112     //
    113     // Distance from source position to center of ellipse.
    114     // If the distance is 0 distance, Alpha is not specified.
    115     // The calculation has failed and returnes kFALSE.
    116     //
    117     const Double_t dist = sqrt(sx*sx + sy*sy);  // [mm]
    118     if (dist==0)
    119         return kFALSE;
    120 
    121     //
    122     // Calculate Alpha and Cosda = cos(d,a)
    123     // The sign of Cosda will be used for quantities containing
    124     // a head-tail information
    125     //
    126     // *OLD* const Double_t arg = (sy-tand*sx) / (dist*sqrt(tand*tand+1));
    127     // *OLD* fAlpha = asin(arg)*kRad2Deg;
    128     //
    129     const Double_t arg1 = cd*sy-sd*sx;          // [mm]
    130     const Double_t arg2 = cd*sx+sd*sy;          // [mm]
    131 
    132     //
    133     // Due to numerical uncertanties in the calculation of the
    134     // square root (dist) and arg1 it can happen (in less than 1e-5 cases)
    135     // that the absolute value of arg exceeds 1. Because this uncertainty
    136     // results in an Delta Alpha which is still less than 1e-3 we don't care
    137     // about this uncertainty in general and simply set values which exceed
    138     // to 1 saving its sign.
    139     //
    140     const Double_t arg = arg1/dist;
    141     fAlpha = TMath::Abs(arg)>1 ? TMath::Sign(90., arg) : asin(arg)*kRad2Deg;  // [deg]
    142 
    143     fCosDeltaAlpha = arg2/dist;                 // [1]
    144     fDist          = dist;                      // [mm]
    145 
    146     SetReadyToSave();
    147 
    148     return kTRUE;
    149 }
     110Int_t MHillasSrc::Calc(const MHillas &hillas)
     111{
     112     const Double_t mx = hillas.GetMeanX();     // [mm]
     113     const Double_t my = hillas.GetMeanY();     // [mm]
     114
     115     const Double_t sx = mx - fSrcPos->GetX();   // [mm]
     116     const Double_t sy = my - fSrcPos->GetY();   // [mm]
     117
     118     const Double_t sd = hillas.GetSinDelta();  // [1]
     119     const Double_t cd = hillas.GetCosDelta();  // [1]
     120
     121     //
     122     // Distance from source position to center of ellipse.
     123     // If the distance is 0 distance, Alpha is not specified.
     124     // The calculation has failed and returnes kFALSE.
     125     //
     126     const Double_t dist = TMath::Sqrt(sx*sx + sy*sy);  // [mm]
     127     if (dist==0)
     128         return 1;
     129
     130     //
     131     // Calculate Alpha and Cosda = cos(d,a)
     132     // The sign of Cosda will be used for quantities containing
     133     // a head-tail information
     134     //
     135     // *OLD* const Double_t arg = (sy-tand*sx) / (dist*sqrt(tand*tand+1));
     136     // *OLD* fAlpha = asin(arg)*kRad2Deg;
     137     //
     138     const Double_t arg2 = cd*sx + sd*sy;          // [mm]
     139     if (arg2==0)
     140         return 2;
     141
     142     const Double_t arg1 = cd*sy - sd*sx;          // [mm]
     143
     144     //
     145     // Due to numerical uncertanties in the calculation of the
     146     // square root (dist) and arg1 it can happen (in less than 1e-5 cases)
     147     // that the absolute value of arg exceeds 1. Because this uncertainty
     148     // results in an Delta Alpha which is still less than 1e-3 we don't care
     149     // about this uncertainty in general and simply set values which exceed
     150     // to 1 saving its sign.
     151     //
     152     const Double_t arg = arg1/dist;
     153     fAlpha = TMath::Abs(arg)>1 ? TMath::Sign(90., arg) : TMath::ASin(arg)*TMath::RadToDeg(); // [deg]
     154
     155     fCosDeltaAlpha = arg2/dist;                 // [1]
     156     fDist          = dist;                      // [mm]
     157
     158     // ----- Calculation of Distance to closest approach 'DCA' -----
     159
     160     // Components of DCA vector
     161     const Double_t fpd1 = sx - arg2*cd;         // [mm]
     162     const Double_t fpd2 = sy - arg2*sd;         // [mm]
     163
     164     // Determine the correct sign of the DCA (cross product of DCA vector and the
     165     // vector going from the intersection point of the DCA vector with the shower axis
     166     // to the COG)
     167     const Double_t sign = arg2*cd*fpd2 - arg2*sd*fpd1;
     168     fDCA  = TMath::Sign(TMath::Sqrt(fpd1*fpd1 + fpd2*fpd2), sign); // [mm]
     169
     170     // Calculate angle of the shower axis with respect to the x-axis
     171     fDCADelta = TMath::ACos((sx-fpd1)/TMath::Abs(arg2))*TMath::RadToDeg(); // [deg]
     172
     173     // Enlarge the interval of fDdca to [0, 2pi]
     174     if (sy < fpd2)
     175         fDCADelta = 360 - fDCADelta;
     176
     177     SetReadyToSave();
     178
     179     return 0;
     180}
    150181
    151182// --------------------------------------------------------------------------
     
    160191    *fLog << " - Alpha          [deg] = " << fAlpha << endl;
    161192    *fLog << " - CosDeltaAlpha        = " << fCosDeltaAlpha << endl;
     193    *fLog << " - DCA            [mm]  = " << fDCA << endl;
     194    *fLog << " - DCA delta      [deg] = " << fDCADelta*TMath::RadToDeg() << endl;
    162195}
    163196
     
    174207    *fLog << " - Alpha          [deg] = " << fAlpha << endl;
    175208    *fLog << " - CosDeltaAlpha        = " << fCosDeltaAlpha << endl;
     209    *fLog << " - DCA            [deg] = " << fDCA*geom.GetConvMm2Deg() << endl;
     210    *fLog << " - DCA delta      [deg] = " << fDCADelta*TMath::RadToDeg() << endl;
    176211}
    177212
     
    183218void MHillasSrc::Set(const TArrayF &arr)
    184219{
    185     if (arr.GetSize() != 3)
     220    if (arr.GetSize() != 5)
    186221        return;
    187222
     
    189224    fDist  = arr.At(1);         // [mm]   distance between src and center of ellipse
    190225    fCosDeltaAlpha = arr.At(2); // [1]    cosine of angle between d and a
    191 }
     226    fDCA = arr.At(3);           // [mm]
     227    fDCADelta = arr.At(4);      // [mm]
     228}
  • trunk/MagicSoft/Mars/mimage/MHillasSrc.h

    r4710 r4817  
    1717    Float_t fCosDeltaAlpha; // [1]    cosine of angle between d and a
    1818
     19    Float_t fDCA;           // [mm]   Distance to closest approach 'DCA'
     20    Float_t fDCADelta;      // [deg]  Angle of the shower axis with respect to the x-axis
     21
    1922public:
    2023    MHillasSrc(const char *name=NULL, const char *title=NULL);
     
    2831    Float_t GetDist()          const { return fDist; }
    2932    Float_t GetCosDeltaAlpha() const { return fCosDeltaAlpha; }
     33    Float_t GetDCA()           const { return fDCA; }
     34    Float_t GetDCADelta()      const { return fDCADelta; }
    3035
    3136    void Print(Option_t *opt=NULL) const;
    3237    void Print(const MGeomCam &geom) const;
    3338
    34     virtual Bool_t Calc(const MHillas &hillas);
     39    virtual Int_t Calc(const MHillas &hillas);
    3540
    3641    void Set(const TArrayF &arr);
    3742
    38     ClassDef(MHillasSrc, 4) // Container to hold source position dependant parameters
     43    ClassDef(MHillasSrc, 5) // Container to hold source position dependant parameters
    3944};
    4045
Note: See TracChangeset for help on using the changeset viewer.