Ignore:
Timestamp:
11/01/04 22:28:46 (20 years ago)
Author:
mase
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars/mmuon
Files:
3 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/mmuon/MMuonCalibPar.h

    r5210 r5331  
    1010#endif
    1111
    12 class MGeomCam;
    1312class MCerPhotEvt;
    1413class MMuonSearchPar;
     
    2928  Float_t fArcPhiThres;   // The threshold value to define arc phi
    3029  Float_t fArcWidthThres; // The threshold value to define arc width
    31   Int_t   fArcPhiBinNum;  // The bin number for the histogram of arc phi. You may change this value. However, if you change this, YOU ALSO HAVE TO CHANGE THE THRESHOLD VALUE TO GET ARC LENGTH.
     30
     31public:
     32  MMuonCalibPar(const char *name=NULL, const char *title=NULL);
     33  ~MMuonCalibPar();
     34 
     35  Int_t   fArcPhiBinNum;         // The bin number for the histogram of arc phi. You may change this value. However, if you change this, YOU ALSO HAVE TO CHANGE THE THRESHOLD VALUE TO GET ARC LENGTH.
    3236  Int_t   fArcWidthBinNum;       // The bin number for the histogram of arc wid
    3337  Float_t fArcPhiHistStartVal;   // The starting value for the histogram of arc phi
     
    3539  Float_t fArcWidthHistStartVal; // The starting value for the histogram of arc width
    3640  Float_t fArcWidthHistEndVal;   // The end value for the histogram of arc width
    37   Bool_t fEnableImpactCalc;      // If true, the impact calculation will be done, which will consume a lot of time.
    38   Bool_t fDisablePreCuts;        // If true, the pre cuts to select muons for the calibration will be disabled.
    39   Bool_t fUseCleanForWidth;      // If true, the only signal after image cleaningwill be used for the histogram of arc width
    4041
    41   TH1F *fHistPhi;      // Histogram of photon distribution along the arc.
    42   TH1F *fHistWidth;    // Histogram of radial photon distribution of the arc.
    43  
    44 public:
    45   MMuonCalibPar(const char *name=NULL, const char *title=NULL);
    46   ~MMuonCalibPar();
     42  TH1F *fHistPhi;   // Histogram of photon distribution along the arc.
     43  TH1F *fHistWidth; // Histogram of radial photon distribution of the arc.
    4744 
    4845  void Reset();
     
    6562  TH1F    *GetHistWidth()     { return fHistWidth; }
    6663 
     64  void    SetArcLength(Float_t len)       { fArcLength = len; }
     65  void    SetArcPhi(Float_t phi)          { fArcPhi = phi; }
     66  void    SetArcWidth(Float_t wid)        { fArcWidth = wid; }
     67  void    SetChiArcPhi(Float_t chi)       { fChiArcPhi = chi; }
     68  void    SetChiArcWidth(Float_t chi)     { fChiArcWidth = chi; }
    6769  void    SetMargin(Float_t margin)       { fMargin = margin; }
     70  void    SetMuonSize(Float_t size)       { fMuonSize = size; }
     71  void    SetEstImpact(Float_t impact)    { fEstImpact = impact; }
     72  void    SetUseUnmap()                   { fUseUnmap = kTRUE; }
     73  void    SetPeakPhi(Float_t phi)         { fPeakPhi = phi; }
    6874  void    SetArcPhiThres(Float_t thres)   { fArcPhiThres = thres; }
    6975  void    SetArcWidthThres(Float_t thres) { fArcWidthThres = thres; }
     
    7177  void    SetArcWidthBinNum(Int_t num)    { fArcWidthBinNum = num; }
    7278 
    73   void    EnableImpactCalc()   { fEnableImpactCalc = kTRUE; }
    74   Bool_t  IsEnableImpactCalc() { return fEnableImpactCalc; }
    75   void    DisablePreCuts()     { fDisablePreCuts = kTRUE; }
    76   Bool_t  IsDisablePreCuts()   { return fDisablePreCuts; }
    77   void    UseCleanForWidth()   { fUseCleanForWidth = kTRUE; }
    78   Bool_t  IsUseCleanForWidth() { return fUseCleanForWidth; }
    79 
    8079  void    Print(Option_t *opt=NULL) const;
    8180
    82   void    FillHist(const MGeomCam &geom, const MCerPhotEvt &evt,
    83                    const MMuonSearchPar &musearch);
    84   void    CalcPhi(const MGeomCam &geom, const MCerPhotEvt &evt,
    85                   const MMuonSearchPar &musearch);
    86   void    CalcImpact(const MGeomCam &geom, const MMuonSearchPar &musearch,
    87                      Int_t effbinnum, Float_t startfitval, Float_t endfitval);
    88   Float_t CalcWidth(const MGeomCam &geom, const MCerPhotEvt &evt,
    89                   const MMuonSearchPar &musearch);
    90   Int_t   Calc(const MGeomCam &geom, const MCerPhotEvt &evt,
    91                MMuonSearchPar &musearch, const Float_t *cuts);
    92  
    9381  ClassDef(MMuonCalibPar, 1) // Container to hold muon calibration parameters
    9482};
  • trunk/MagicSoft/Mars/mmuon/MMuonSearchPar.cc

    r5210 r5331  
    3030// Storage Container for muon
    3131//
    32 //  This class is the container for muon parameters. Actually, in this class,
    33 // muons are searched by fitting the image with a circle. (This function will
    34 // be called by using the class of MMuonSearchParCalc.) This container
    35 // especially holds the information of the results of the search (fit of
    36 // a image by a circle).
     32//  This class is the container for muon parameters. Actually, the calculation
     33// is done here. Muons are searched by fitting the image with a circle.
     34// (This function will be called by using the class of MMuonSearchParCalc.)
     35// This container especially holds the information of the results of the
     36// search (fit of a image by a circle).
    3737//
    3838//  In order to use further information of muons such as the width of arcs,
     
    8585#include "MCerPhotEvt.h"
    8686#include "MCerPhotPix.h"
    87 #include "MBinning.h"
    8887
    8988using namespace std;
     
    118117//
    119118void MMuonSearchPar::CalcTempCenter(const MHillas &hillas,
    120       Float_t *xtmp1, Float_t *ytmp1, Float_t *xtmp2, Float_t *ytmp2)
     119      Float_t &xtmp1, Float_t &ytmp1, Float_t &xtmp2, Float_t &ytmp2)
    121120{
    122121  Float_t a,dx,dy;
     
    128127  dy = -tmp_r/TMath::Sqrt(1+a*a)/3.;
    129128
    130   *xtmp1 = hillas.GetMeanX() + dx;
    131   *ytmp1 = hillas.GetMeanY() + dy;
    132   *xtmp2 = hillas.GetMeanX() - dx;
    133   *ytmp2 = hillas.GetMeanY() - dy;
     129  xtmp1 = hillas.GetMeanX() + dx;
     130  ytmp1 = hillas.GetMeanY() + dy;
     131  xtmp2 = hillas.GetMeanX() - dx;
     132  ytmp2 = hillas.GetMeanY() - dy;
    134133}
    135134
     
    140139//
    141140Bool_t MMuonSearchPar::CalcRadius(const MGeomCam &geom, const MCerPhotEvt &evt,
    142       Float_t x, Float_t y, Float_t *r, Float_t *sigma)
     141      Float_t x, Float_t y, Float_t &r, Float_t &sigma)
    143142{
    144143  Float_t mean_r=0., dev_r=0., sums=0., tmp=0.;
     
    165164      return kFALSE;
    166165
    167   *r = mean_r/sums;
    168 
    169   if(dev_r/sums-(*r)*(*r)<1.E-10)
     166  r = mean_r/sums;
     167
     168  if(dev_r/sums-(r)*(r)<1.E-10)
    170169      return kFALSE;
    171170
    172   *sigma = TMath::Sqrt(dev_r/sums-(*r)*(*r));
     171  sigma = TMath::Sqrt(dev_r/sums-(r)*(r));
    173172
    174173  return kTRUE;
     
    182181void MMuonSearchPar::CalcMinimumDeviation
    183182( const MGeomCam &geom, const MCerPhotEvt &evt, Float_t x, Float_t y,
    184  Float_t xcog, Float_t ycog, Float_t sigma, Float_t *opt_rad,
    185  Float_t *new_sigma, Float_t *newx, Float_t *newy )
     183 Float_t xcog, Float_t ycog, Float_t sigma, Float_t &opt_rad,
     184 Float_t &new_sigma, Float_t &newx, Float_t &newy )
    186185{
    187186  Float_t delta = 3.;  // 3 mm (1/10 of an inner pixel size) Step to move.
     
    195194      if(r2 > 360000.)
    196195      {
    197           *new_sigma=sigma;
    198           *opt_rad=rad_tmp;
    199           *newx=x; *newy=y;
     196          new_sigma=sigma;
     197          opt_rad=rad_tmp;
     198          newx=x;
     199          newy=y;
    200200          fNoMuon = kTRUE;
    201201          break;
    202202      }
    203       if(CalcRadius(geom,evt,x,y+delta,&rad_tmp,&sig_tmp))
     203      if(CalcRadius(geom,evt,x,y+delta,rad_tmp,sig_tmp))
    204204      {
    205205          if(sig_tmp<sigma)
    206206          {
    207               sigma=sig_tmp; *opt_rad=rad_tmp;
     207              sigma=sig_tmp;
     208              opt_rad=rad_tmp;
    208209              y=y+delta;
    209210              continue;
    210211          }
    211212      }
    212       if(CalcRadius(geom,evt,x-delta,y,&rad_tmp,&sig_tmp))
     213      if(CalcRadius(geom,evt,x-delta,y,rad_tmp,sig_tmp))
    213214      {
    214215          if(sig_tmp<sigma)
    215216          {
    216               sigma=sig_tmp; *opt_rad=rad_tmp;
     217              sigma=sig_tmp;
     218              opt_rad=rad_tmp;
    217219              x=x-delta;
    218220              continue;
    219221          }
    220222      }
    221       if(CalcRadius(geom,evt,x+delta,y,&rad_tmp,&sig_tmp))
     223      if(CalcRadius(geom,evt,x+delta,y,rad_tmp,sig_tmp))
    222224      {
    223225          if(sig_tmp<sigma)
    224226          {
    225               sigma=sig_tmp; *opt_rad=rad_tmp;
     227              sigma=sig_tmp;
     228              opt_rad=rad_tmp;
    226229              x=x+delta;
    227230              continue;
    228231          }
    229232      }
    230       if(CalcRadius(geom,evt,x,y-delta,&rad_tmp,&sig_tmp))
     233      if(CalcRadius(geom,evt,x,y-delta,rad_tmp,sig_tmp))
    231234      {
    232235          if(sig_tmp<sigma)
    233236          {
    234               sigma=sig_tmp; *opt_rad=rad_tmp;
     237              sigma=sig_tmp;
     238              opt_rad=rad_tmp;
    235239              y=y-delta;
    236240              continue;
    237241          }
    238242      }
    239       *new_sigma=sigma;
    240       *newx=x; *newy=y;
     243      new_sigma=sigma;
     244      newx=x;
     245      newy=y;
    241246      break;
    242247    }
     
    257262   
    258263  // gets temporaly center
    259   CalcTempCenter(hillas,&xtmp1,&ytmp1,&xtmp2,&ytmp2);
     264  CalcTempCenter(hillas,xtmp1,ytmp1,xtmp2,ytmp2);
    260265 
    261266  // determine which position will be the true position. Here mainly
    262267  // the curvature of a muon arc is relied on.
    263   CalcRadius(geom, evt, xtmp1,ytmp1,&rad,&dev);
    264   CalcRadius(geom, evt, xtmp2,ytmp2,&rad2,&dev2);
     268  CalcRadius(geom, evt, xtmp1,ytmp1,rad,dev);
     269  CalcRadius(geom, evt, xtmp2,ytmp2,rad2,dev2);
    265270  if(dev2<dev){
    266271    xtmp1=xtmp2; ytmp1=ytmp2; dev=dev2; rad=rad2;
     
    268273 
    269274  // find the best fit.
    270   CalcMinimumDeviation(geom, evt, xtmp1,ytmp1,hillas.GetMeanX(),hillas.GetMeanY(),
    271                  dev, &opt_rad, &new_sigma, &newx, &newy);
     275  CalcMinimumDeviation(geom, evt, xtmp1,ytmp1,hillas.GetMeanX(),
     276                       hillas.GetMeanY(), dev, opt_rad, new_sigma,
     277                       newx, newy);
    272278
    273279  fRadius = opt_rad;
  • trunk/MagicSoft/Mars/mmuon/MMuonSearchPar.h

    r5210 r5331  
    44#ifndef MARS_MParContainer
    55#include "MParContainer.h"
    6 #endif
    7 
    8 #ifndef ROOT_TH1
    9 #include <TH1.h>
    10 #endif
    11 
    12 #ifndef ROOT_TArrayI
    13 #include <TArrayI.h>
    146#endif
    157
     
    3931    void    SetNoMuon()          { fNoMuon = kTRUE; }
    4032
     33    void   CalcTempCenter(const MHillas &hillas, Float_t &xtmp1,
     34                          Float_t &ytmp1, Float_t &xtmp2, Float_t &ytmp2);
     35    Bool_t CalcRadius(const MGeomCam &geom, const MCerPhotEvt &evt, Float_t x,
     36                      Float_t y, Float_t &r, Float_t &sigma);
     37    void   CalcMinimumDeviation(const MGeomCam &geom, const MCerPhotEvt &evt,
     38                             Float_t x, Float_t y, Float_t xcog,
     39                             Float_t ycog, Float_t sigma, Float_t &opt_rad,
     40                             Float_t &new_sigma, Float_t &newx, Float_t &newy);
     41    void   Calc(const MGeomCam &geom, const MCerPhotEvt &evt,
     42                const MHillas &hillas);
     43
    4144    void   Print(Option_t *opt=NULL) const;
    4245    void   Print(const MGeomCam &geom, Option_t *opt=NULL) const;
    43     void   CalcTempCenter(const MHillas &hillas, Float_t *xtmp1,
    44                          Float_t *ytmp1, Float_t *xtmp2, Float_t *ytmp2);
    45     Bool_t CalcRadius(const MGeomCam &geom, const MCerPhotEvt &evt, Float_t x,
    46                       Float_t y, Float_t *r, Float_t *sigma);
    47     void   CalcMinimumDeviation(const MGeomCam &geom, const MCerPhotEvt &evt,
    48                          Float_t x, Float_t y, Float_t xcog,
    49                          Float_t ycog, Float_t sigma, Float_t *opt_rad,
    50                          Float_t *new_sigma, Float_t *newx, Float_t *newy);
    51     void   Calc(const MGeomCam &geom, const MCerPhotEvt &evt,
    52                const MHillas &hillas);
    5346
    5447    ClassDef(MMuonSearchPar, 1) // Container to hold muon search parameters
Note: See TracChangeset for help on using the changeset viewer.