Changeset 7365


Ignore:
Timestamp:
09/26/05 19:17:18 (19 years ago)
Author:
meyer
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • trunk/MagicSoft/Mars/Changelog

    r7363 r7365  
    1818
    1919                                                 -*-*- END OF LINE -*-*-
     20
     21 2005/09/26 Markus Meyer
     22
     23   * mmuon/MMuonSearchPar.cc/.h:
     24     - Mean arrival time of core Pixels is calculated
     25
     26   * mmuon/MHSingleMuon.cc/.h:
     27     - new histogram fHistTime (TH1), which contents difference
     28       of the arrival time between Mean arrival time of the core
     29       pixels and every pixel in the margin around the circle.
     30       The histogram is fitted by an gaussian distribution.
     31     - fHistPhi is now a TProfile
     32     - fHistPhi is only filled with pixels with an arrival time 
     33       of < 2 sigma difference to the fitted mean value
     34
     35   * mmuon/MMuonSetup.cc:
     36     - fThresholdArcPhi is changed to 5 phe
     37   
     38   * mmuon/MHMuonPar.cc:
     39     - for the reference values to MC, the integral over the mean
     40       values of the different radius bins is taken again
     41
     42
    2043
    2144 2005/09/22 Daniela Dorner
  • trunk/MagicSoft/Mars/mmuon/MHMuonPar.cc

    r7253 r7365  
    258258    const Int_t bin2 = p.GetXaxis()->FindFixBin(b);
    259259
    260     return Integral(p, bin1, bin2);
     260//    return Integral(p, bin1, bin2);
     261    return p.Integral(bin1, bin2);
    261262}
    262263
     
    265266    if (TString(opt)==TString("pad4"))
    266267    {
    267         const TString txt = Form("\\Sigma_{%.2f\\circ}^{%.2f\\circ} = %3.1fm",
    268                                  fgIntegralLoLim, fgIntegralUpLim, Integral(fHistBroad)*1000);
     268        const TString txt = Form("\\Sigma_{%.2f\\circ}^{%.2f\\circ} = %.3f",
     269                                 fgIntegralLoLim, fgIntegralUpLim, Integral(fHistBroad));
     270//                                 fgIntegralLoLim, fgIntegralUpLim, Integral(fHistBroad)*1000);
    269271
    270272        TLatex text(0.55, 0.93, txt);
  • trunk/MagicSoft/Mars/mmuon/MHSingleMuon.cc

    r7225 r7365  
    127127    fHistWidth.UseCurrentStyle();
    128128
     129    fHistTime.SetName("HistTime");
     130    fHistTime.SetTitle("HistTime");
     131    fHistTime.SetXTitle("timing difference");
     132    fHistTime.SetYTitle("number of pixels");
     133    fHistTime.SetDirectory(NULL);
     134    fHistTime.SetFillStyle(4000);
     135    fHistTime.UseCurrentStyle();
     136
    129137    MBinning bins;
    130138    bins.SetEdges(20, -180, 180);
     
    133141    bins.SetEdges(28, 0.3, 1.7);
    134142    bins.Apply(fHistWidth);
     143
     144    bins.SetEdges(100, -10, 10);
     145    bins.Apply(fHistTime);
    135146}
    136147
     
    169180    ApplyBinning(*plist, "ArcPhi",    &fHistPhi);
    170181    ApplyBinning(*plist, "MuonWidth", &fHistWidth);
     182    ApplyBinning(*plist, "MuonTime", &fHistTime);
    171183
    172184    return kTRUE;
     
    182194    fHistPhi.Reset();
    183195    fHistWidth.Reset();
     196    fHistTime.Reset();
    184197
    185198    const Int_t entries = fSignalCam->GetNumPixels();
     
    202215        if (dist < fMuonSearchPar->GetRadius() + fMargin &&
    203216            dist > fMuonSearchPar->GetRadius() - fMargin)
    204             fHistPhi.Fill(TMath::ATan2(dx, dy)*TMath::RadToDeg(), pix.GetNumPhotons());
    205 
     217        {
     218            fHistTime.Fill(pix.GetArrivalTime()-fMuonSearchPar->GetTime());
     219        }
    206220        // use only the inner pixles. FIXME: This is geometry dependent
    207221        if(i>397)
     
    210224        fHistWidth.Fill(dist*fGeomCam->GetConvMm2Deg(), pix.GetNumPhotons());
    211225    }
    212 
     226    // Setup the function and perform the fit
     227    TF1 g1("g1", "gaus", -10, 10);
     228
     229    // Choose starting values as accurate as possible
     230    g1.SetParameter(0, fHistTime.GetMaximum());
     231    g1.SetParameter(1, 0);
     232    g1.SetParameter(2, 0.3);
     233    g1.SetParLimits(1, -0.5, 0.5);
     234    g1.SetParLimits(2, 0, 1);
     235    // options : N  do not store the function, do not draw
     236    //           I  use integral of function in bin rather than value at bin center
     237    //           R  use the range specified in the function range
     238    //           Q  quiet mode
     239    fHistTime.Fit(&g1, "QNRB");
     240
     241    //    Double_t err;
     242    Double_t mean, meanerr, sig, sigerr;
     243    gMinuit->GetParameter(2, sig, sigerr); // get the sigma value
     244    gMinuit->GetParameter(1, mean, meanerr); // get the sigma value
     245
     246    for (Int_t i=0; i<entries; i++)
     247    {
     248        const MSignalPix &pix  = (*fSignalCam)[i];
     249        const MGeomPix   &gpix = (*fGeomCam)[i];
     250
     251        const Float_t dx = gpix.GetX() - cenx;
     252        const Float_t dy = gpix.GetY() - ceny;
     253
     254        const Float_t dist = TMath::Hypot(dx, dy);
     255
     256        // if the signal is not near the estimated circle, it is ignored.
     257        if (dist < fMuonSearchPar->GetRadius() + fMargin &&
     258            dist > fMuonSearchPar->GetRadius() - fMargin)
     259        {
     260            if(TMath::Abs(pix.GetArrivalTime()-(fMuonSearchPar->GetTime()+mean))<2*sig)
     261                fHistPhi.Fill(TMath::ATan2(dx, dy)*TMath::RadToDeg(), pix.GetNumPhotons());
     262        }
     263    }
    213264    // error estimation (temporaly)
    214265    //  The error is estimated from the signal. In order to do so, we have to
     
    244295// below the threshold is found, kFALSE is returned.
    245296//
    246 Bool_t MHSingleMuon::FindRangeAboveThreshold(const TH1 &h, Float_t thres, Int_t &first, Int_t &last) const
     297Bool_t MHSingleMuon::FindRangeAboveThreshold(const TProfile &h, Float_t thres, Int_t &first, Int_t &last) const
    247298{
    248299    const Int_t n      = h.GetNbinsX();
     
    362413    //           R  use the range specified in the function range
    363414    //           Q  quiet mode
    364     fHistWidth.Fit(&f1, "QR0");
     415//    fHistWidth.Fit(&f1, "QRO");
     416    fHistWidth.Fit(&f1, "QRN");
    365417
    366418    chi = f1.GetChisquare()/f1.GetNDF();
     
    477529    gPad->SetBorderMode(0);
    478530    fHistWidth.Draw();
    479 
    480 }
     531}
  • trunk/MagicSoft/Mars/mmuon/MHSingleMuon.h

    r7225 r7365  
    2424    Double_t fMargin;               //!
    2525
    26     TH1F    fHistPhi;   // Histogram of photon distribution along the arc.
     26    TProfile fHistPhi;   // Histogram of photon distribution along the arc.
    2727    TProfile fHistWidth; // Histogram of radial photon distribution of the arc.
     28    TH1F     fHistTime;   // Histogram of arrival time distribution along the arc.
    2829
    29     Bool_t FindRangeAboveThreshold(const TH1 &h, Float_t thres, Int_t &first, Int_t &last) const;
     30    Bool_t FindRangeAboveThreshold(const TProfile &h, Float_t thres, Int_t &first, Int_t &last) const;
    3031
    3132public:
     
    3839    Bool_t CalcWidth(Double_t, Double_t &, Double_t &);
    3940
    40     const TH1F &GetHistPhi() const { return fHistPhi; }
     41    const TProfile &GetHistPhi() const { return fHistPhi; }
     42    const TH1F     &GetHistTime() const { return fHistTime; }
    4143    const TProfile &GetHistWidth() const { return fHistWidth; }
    4244
  • trunk/MagicSoft/Mars/mmuon/MMuonSearchPar.cc

    r7176 r7365  
    8989    fCenterX   =  0;
    9090    fCenterY   =  0;
     91    fTime      =  0;
    9192}
    9293
     
    150151
    151152    Int_t p=0;
     153    Int_t q=0;
     154
     155    Float_t fTime2=0;
     156
    152157    for (int i=0; i<n; i++)
    153158    {
    154159        const MSignalPix &pix = evt[i];
     160
    155161        if (pix.IsPixelUsed())
    156162        {
     
    160166            fY[p] = geom[i].GetY();
    161167            p++;
     168
     169            //timing
     170            if(pix.IsPixelCore())
     171            {
     172                fTime  += pix.GetArrivalTime();
     173                fTime2 += pix.GetArrivalTime()*pix.GetArrivalTime();
     174                q++;
     175            }
    162176        }
    163177    }
     178    fTime    = fTime/q;
     179    fTime2   = fTime2/q;
     180    fTimeRms = TMath::Sqrt(fTime2-(fTime*fTime));
     181
    164182    fSignal.Set(p);
    165183
  • trunk/MagicSoft/Mars/mmuon/MMuonSearchPar.h

    r7134 r7365  
    2121    Float_t fCenterX;   // An estimated center position in X of the muon ring [mm]
    2222    Float_t fCenterY;   // An estimated center position in Y of the muon ring [mm]
     23    Float_t fTime;
     24    Float_t fTimeRms;
    2325
    2426    MArrayF fSignal;    //! Temporary storage for signal
     
    4042    Float_t GetCenterX()   const { return fCenterX; }
    4143    Float_t GetCenterY()   const { return fCenterY; }
     44    Float_t GetTime()      const { return fTime; }
     45    Float_t GetTimeRms()   const { return fTimeRms; }
    4246
    4347    // MMuonSearchPar
  • trunk/MagicSoft/Mars/mmuon/MMuonSetup.cc

    r7225 r7365  
    4949//
    5050MMuonSetup::MMuonSetup(const char *name, const char *title)
    51     : fMargin(0.2), fThresholdArcPhi(30), fThresholdArcWidth(2)
     51    : fMargin(0.2), fThresholdArcPhi(5), fThresholdArcWidth(2)
    5252{
    5353    fName  = name  ? name  : "MMuonSetup";
Note: See TracChangeset for help on using the changeset viewer.