Changeset 7538 for trunk/MagicSoft


Ignore:
Timestamp:
02/27/06 13:04:28 (19 years ago)
Author:
tbretz
Message:
*** empty log message ***
Location:
trunk/MagicSoft/Mars
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • trunk/MagicSoft/Mars/Changelog

    r7537 r7538  
    4949     - increased class version of MNewImagePar
    5050
     51   * mimage/MHillasExt.cc:
     52     - added some comments
     53     - removed an obsolete TMath::Abs around "dist"
     54
    5155   * mranforest/MRanForest.[h,cc]:
    5256     - the initialization of fTreeHad was done at the wrong moment

  • trunk/MagicSoft/Mars/mimage/MHNewImagePar.cc

    r6977 r7538  
    1919!   Author(s): Thomas Bretz, 04/2003 <mailto:tbretz@astro.uni-wuerzburg.de>
    2020!
    21 !   Copyright: MAGIC Software Development, 2000-2004
     21!   Copyright: MAGIC Software Development, 2000-2006
    2222!
    2323!
    … …  
    2727//
    2828// MHNewImagePar
     29//
     30//
     31// Version 2:
     32// ----------
     33//  + TH1F fHistConcCOG;   // [ratio] concentration around the center of gravity (all pixels)
     34//  + TH1F fHistConcCore;  // [ratio] concentration of signals inside ellipse (used pixels)
    2935//
    3036////////////////////////////////////////////////////////////////////////////
    … …  
    116122    fHistConc.SetDirectory(NULL);
    117123    fHistConc1.SetDirectory(NULL);
     124    fHistConcCOG.SetDirectory(NULL);
     125    fHistConcCore.SetDirectory(NULL);
    118126    fHistConc.SetName("Conc2");
    119127    fHistConc1.SetName("Conc1");
     128    fHistConcCOG.SetName("ConcCOG");
     129    fHistConcCore.SetName("ConcCore");
    120130    fHistConc.SetTitle("Ratio: Conc");
    121131    fHistConc1.SetTitle("Ratio: Conc1");
     132    fHistConcCOG.SetTitle("Ratio: ConcCOG");
     133    fHistConcCore.SetTitle("Ratio: ConcCore");
    122134    fHistConc.SetXTitle("Ratio");
    123135    fHistConc1.SetXTitle("Ratio");
     136    fHistConcCOG.SetXTitle("Ratio");
     137    fHistConcCore.SetXTitle("Ratio");
    124138    fHistConc.SetYTitle("Counts");
    125139    fHistConc1.SetYTitle("Counts");
     140    fHistConcCOG.SetYTitle("Counts");
     141    fHistConcCore.SetYTitle("Counts");
     142    fHistConc.UseCurrentStyle();
    126143    fHistConc1.UseCurrentStyle();
    127     fHistConc.UseCurrentStyle();
     144    fHistConcCOG.UseCurrentStyle();
     145    fHistConcCore.UseCurrentStyle();
    128146    fHistConc.SetFillStyle(4000);
    129147    fHistConc1.SetFillStyle(4000);
     148    fHistConcCOG.SetFillStyle(4000);
     149    fHistConcCore.SetFillStyle(4000);
    130150    fHistConc1.SetLineColor(kBlue);
    131 
     151    fHistConcCOG.SetLineColor(kBlue);
    132152
    133153    MBinning bins;
    … …  
    138158    bins.Apply(fHistConc);
    139159    bins.Apply(fHistConc1);
     160    bins.Apply(fHistConcCOG);
     161    bins.Apply(fHistConcCore);
    140162
    141163    bins.SetEdges(75, 0.5, 150.5);
    … …  
    189211    //ApplyBinning(*plist, "Area",    &fHistCoreArea);
    190212
    191     ApplyBinning(*plist, "Conc",    &fHistConc);
    192     ApplyBinning(*plist, "Conc1",   &fHistConc1);
     213    ApplyBinning(*plist, "Conc",     &fHistConc);
     214    ApplyBinning(*plist, "Conc1",    &fHistConc1);
     215    ApplyBinning(*plist, "ConcCOG",  &fHistConcCOG);
     216    ApplyBinning(*plist, "ConcCore", &fHistConcCore);
    193217
    194218    return kTRUE;
    … …  
    223247    fHistConc.Fill(h.GetConc(), w);
    224248    fHistConc1.Fill(h.GetConc1(), w);
     249    fHistConcCOG.Fill(h.GetConcCOG(), w);
     250    fHistConcCore.Fill(h.GetConcCore(), w);
    225251
    226252    return kTRUE;
    … …  
    284310void MHNewImagePar::Paint(Option_t *o)
    285311{
    286     if (fHistLeakage1.GetMaximum()>0 && gPad->GetPad(1))
    287         gPad->GetPad(1)->SetLogy();
     312    if (fHistLeakage1.GetMaximum()>0 && gPad->GetPad(1) && gPad->GetPad(1)->GetPad(1))
     313        gPad->GetPad(1)->GetPad(1)->SetLogy();
    288314}
    289315
    … …  
    306332
    307333    if (!same)
    308         pad->Divide(2,2);
     334        pad->Divide(2, 1);
    309335    else
    310336    {
    … …  
    315341        fHistUsedArea.SetName("UsedAreaSame");
    316342        fHistCoreArea.SetName("CoreAreaSame");
     343        fHistConcCOG.SetName("ConcCOGSame");
     344        fHistConcCore.SetName("ConcCoreSame");
    317345        fHistConc1.SetName("Conc1Same");
    318346        fHistConc.SetName("Conc2Same");
    … …  
    324352        fHistUsedArea.SetDirectory(0);
    325353        fHistCoreArea.SetDirectory(0);
     354        fHistConcCOG.SetDirectory(0);
     355        fHistConcCore.SetDirectory(0);
    326356        fHistConc1.SetDirectory(0);
    327357        fHistConc.SetDirectory(0);
    … …  
    332362        fHistUsedPix.SetLineColor(kCyan);
    333363        fHistConc1.SetLineColor(kMagenta);
     364        fHistConcCOG.SetLineColor(kCyan);
     365        fHistConcCore.SetLineColor(kCyan);
    334366        fHistConc.SetLineColor(kCyan);
    335367        fHistCoreArea.SetLineColor(kMagenta);
    … …  
    338370
    339371    pad->cd(1);
     372    TVirtualPad *pad1=gPad;
     373    pad1->SetBorderMode(0);
     374    pad1->Divide(1,3, 0.001, 0.001);
     375
     376    pad1->cd(1);
    340377    gPad->SetBorderMode(0);
     378    gPad->SetGridx();
    341379    TAxis &x = *fHistLeakage1.GetXaxis();
    342380    x.SetRangeUser(0.0, x.GetXmax());
    … …  
    348386    fHistLeakage2.SetMaximum(0.1);   // dummy value to allow log-scale
    349387
    350     pad->cd(2);
     388    pad1->cd(2);
    351389    gPad->SetBorderMode(0);
     390    gPad->SetGridx();
    352391    RemoveFromPad("UsedPixSame");
    353392    RemoveFromPad("CorePixSame");
    354393    MH::DrawSame(fHistCorePix, fHistUsedPix, "Number of core/used Pixels", same);
    355394
    356     pad->cd(3);
     395    pad1->cd(3);
    357396    gPad->SetBorderMode(0);
     397    gPad->SetGridx();
     398    RemoveFromPad("CoreAreaSame");
     399    RemoveFromPad("UsedAreaSame");
     400    MH::DrawSame(fHistCoreArea, fHistUsedArea, "Area of core/used Pixels", same);
     401
     402    pad->cd(2);
     403    TVirtualPad *pad2=gPad;
     404    pad2->SetBorderMode(0);
     405    pad2->Divide(1, 2, 0.001, 0.001);
     406
     407    pad2->cd(1);
     408    gPad->SetBorderMode(0);
     409    gPad->SetGridx();
    358410    RemoveFromPad("Conc1Same");
    359411    RemoveFromPad("Conc2Same");
    360412    MH::DrawSame(fHistConc1, fHistConc, "Concentrations", same);
    361413
    362     pad->cd(4);
     414    pad2->cd(2);
    363415    gPad->SetBorderMode(0);
    364     RemoveFromPad("CoreAreaSame");
    365     RemoveFromPad("UsedAreaSame");
    366     MH::DrawSame(fHistCoreArea, fHistUsedArea, "Area of core/used Pixels", same);
     416    gPad->SetGridx();
     417    RemoveFromPad("ConcCOGSame");
     418    RemoveFromPad("ConcCoreSame");
     419    MH::DrawSame(fHistConcCore, fHistConcCOG, "Concentrations", same);
    367420}
    368421
    … …  
    373426    if (name.Contains("Leakage2", TString::kIgnoreCase))
    374427        return const_cast<TH1F*>(&fHistLeakage2);
     428    if (name.Contains("ConcCOG", TString::kIgnoreCase))
     429        return const_cast<TH1F*>(&fHistConcCOG);
     430    if (name.Contains("ConcCore", TString::kIgnoreCase))
     431        return const_cast<TH1F*>(&fHistConcCore);
    375432    if (name.Contains("Conc1", TString::kIgnoreCase)) // must be first!
    376433        return const_cast<TH1F*>(&fHistConc1);

  • trunk/MagicSoft/Mars/mimage/MHNewImagePar.h

    r6977 r7538  
    2525    TH1F fHistConc;      // [ratio] concentration ratio: sum of the two highest pixels / fSize
    2626    TH1F fHistConc1;     // [ratio] concentration ratio: sum of the highest pixel / fSize
     27    TH1F fHistConcCOG;   // [ratio] concentration of the three pixels next to COG
     28    TH1F fHistConcCore;  // [ratio] concentration of signals inside or touching the ellipse
    2729
    2830    Float_t fMm2Deg;
    … …  
    5355    TH1F &GetHistConc()      { return fHistConc; }
    5456    TH1F &GetHistConc1()     { return fHistConc1; }
     57    TH1F &GetHistConcCOG()   { return fHistConcCOG; }
     58    TH1F &GetHistConcCore()  { return fHistConcCore; }
    5559
    5660    void SetMmScale(Bool_t mmscale=kTRUE);
    … …  
    6064    void Paint(Option_t *opt="");
    6165
    62     ClassDef(MHNewImagePar, 1) // Histograms of new image parameters
     66    ClassDef(MHNewImagePar, 2) // Histograms of new image parameters
    6367};
    6468

  • trunk/MagicSoft/Mars/mimage/MHillasExt.cc

    r7181 r7538  
    161161
    162162        const Double_t dist = dx*dx+dy*dy;
    163         if (TMath::Abs(dist)>TMath::Abs(maxdist))
     163        if (dist>TMath::Abs(maxdist))
    164164            maxdist = dzx<0 ? -dist : dist;                    // [mm^2]
    165165
    … …  
    177177            maxpix   = nphot;                                  // [1]
    178178            maxpixid = i;//pixid;
    179             continue;                                          // [1]
    180179        }
    181180    }
    … …  
    183182    const MGeomPix &maxp = geom[maxpixid];
    184183
     184    //
     185    // Asymmetry
     186    //
    185187    fAsym = (hil.GetMeanX()-maxp.GetX())*hil.GetCosDelta() +
    186188            (hil.GetMeanY()-maxp.GetY())*hil.GetSinDelta();            // [mm]
    … …  
    195197    fM3Trans = m3y<0 ? -pow(-m3y, 1./3) : pow(m3y, 1./3);      // [mm]
    196198
     199    //
     200    // Distance between max signal and COG
     201    //
    197202    const Double_t md = TMath::Sqrt(TMath::Abs(maxdist));
    198203    fMaxDist = maxdist<0 ? -md : md;                           // [mm]

  • trunk/MagicSoft/Mars/mimage/MNewImagePar.cc

    r7188 r7538  
    3838//    Float_t fConc;                 // [ratio] concentration ratio: sum of the two highest pixels / fSize
    3939//    Float_t fConc1;                // [ratio] concentration ratio: sum of the highest pixel / fSize
     40//    Float_t fConcCOG;              // [ratio] concentration of the three pixels next to COG
     41//    Float_t fConcCore;             // [ratio] concentration of signals inside or touching the ellipse
    4042//
    4143//    Float_t fUsedArea;             // Area of pixels which survived the image cleaning
    … …  
    6567//    + removed fNumSaturatedPixels
    6668//
     69// Version 5:
     70// ----------
     71//  - added fConcCOG
     72//  - added fConcCore
     73//
    6774//
    6875/////////////////////////////////////////////////////////////////////////////
    … …  
    105112    fInnerLeakage1 = -1;
    106113    fInnerLeakage2 = -1;
    107     fInnerSize = -1;
    108 
    109     fConc  = -1;
    110     fConc1 = -1;
     114    fInnerSize     = -1;
     115
     116    fConc     = -1;
     117    fConc1    = -1;
     118    fConcCOG  = -1;
     119    fConcCore = -1;
    111120
    112121    fNumUsedPixels = -1;
    … …  
    131140
    132141    fInnerSize = 0;
     142    fConcCore  = 0;
    133143
    134144    Double_t edgepix1 = 0;
    … …  
    140150    Float_t maxpix1 = 0;                                 // [#phot]
    141151    Float_t maxpix2 = 0;                                 // [#phot]
     152
     153    const Double_t d = geom.GetMaxRadius()*geom.GetMaxRadius();
     154    Double_t dist[3] = { d, d, d };
     155    Int_t    idx[3]  = { -1, -1, -1};
     156
     157    const Double_t rl = 1./(hillas.GetLength()*hillas.GetLength());
     158    const Double_t rw = 1./(hillas.GetWidth() *hillas.GetWidth());
    142159
    143160    const Bool_t ismagiclike =
    … …  
    149166    for (UInt_t i=0; i<npix; i++)
    150167    {
     168        // Get geometry of pixel
     169        const MGeomPix &gpix = geom[i];
     170
     171        // Find the three pixels which are next to the COG
     172        const Double_t dx = gpix.GetX() - hillas.GetMeanX();      // [mm]
     173        const Double_t dy = gpix.GetY() - hillas.GetMeanY();      // [mm]
     174
     175        const Double_t dist0 = dx*dx+dy*dy;
     176
     177        if (dist0<dist[0])
     178        {
     179            dist[2] = dist[1];
     180            dist[1] = dist[0];
     181            dist[0] = dist0;
     182
     183            idx[2]  = idx[1];
     184            idx[1]  = idx[0];
     185            idx[0]  = i;
     186        }
     187        else
     188            if (dist0<dist[1])
     189            {
     190                dist[2] = dist[1];
     191                dist[1] = dist0;
     192
     193                idx[2]  = idx[1];
     194                idx[1]  = i;
     195            }
     196            else
     197                if (dist0<dist[2])
     198                {
     199                    dist[2] = dist0;
     200                    idx[2]  = i;
     201                }
     202
    151203        const MSignalPix &pix = evt[i];
    152204        if (!pix.IsPixelUsed())
    … …  
    156208        if (island>=0 && pix.GetIdxIsland()!=island)
    157209            continue;
    158 
    159         // Get geometry of pixel
    160         //const Int_t pixid = pix->GetPixId();
    161         const MGeomPix &gpix = geom[i/*pixid*/];
    162210
    163211        // count used and core pixels
    … …  
    172220        fUsedArea += gpix.GetA();
    173221
     222        // signal in pixel
    174223        Double_t nphot = pix.GetNumPhotons();
     224
     225        //
     226        // Calculate signal contained inside ellipse
     227        //
     228        const Double_t dzx   =  hillas.GetCosDelta()*dx + hillas.GetSinDelta()*dy; // [mm]
     229        const Double_t dzy   = -hillas.GetSinDelta()*dx + hillas.GetCosDelta()*dy; // [mm]
     230        const Double_t dz    =  gpix.GetD()*gpix.GetD()/4;
     231        const Double_t tana  =  dzy*dzy/(dzx*dzx);
     232        const Double_t distr =  (1+tana)/(rl + tana*rw);
     233        if (distr>dist0-dz || dzx==0)
     234            fConcCore += nphot;
    175235
    176236        //
    … …  
    213273            maxpix2  = maxpix1;
    214274            maxpix1  = nphot;                            // [1]
    215             continue;                                    // [1]
    216275        }
    217 
    218         if (nphot>maxpix2)
    219             maxpix2 = nphot;                             // [1]
     276        else
     277            if (nphot>maxpix2)
     278                maxpix2 = nphot;                         // [1]
    220279    }
    221    
     280
    222281    fInnerLeakage1 = edgepixin1 / fInnerSize;
    223282    fInnerLeakage2 = edgepixin2 / fInnerSize;
    … …  
    226285    fLeakage2 = edgepix2 / hillas.GetSize();
    227286
    228     // FIXME?: in case the pixel with highest signal density is an outer pixel,
    229     // the value of fConc (ratio of signal in two highest pixels to SIZE) should
    230     // rather be 2*fConc1, under the simplest assumption that the light density
    231     // inside the outer (large) pixel is uniform.
     287    // FIXME?: in case the pixel with highest signal density is an outer
     288    // pixel, the value of fConc (ratio of signal in two highest pixels
     289    // to SIZE) should rather be 2*fConc1, under the simplest assumption
     290    // that the light density inside the outer (large) pixel is uniform.
    232291    fConc  = (maxpix1+maxpix2)/hillas.GetSize();         // [ratio]
    233292    fConc1 = maxpix1/hillas.GetSize();                   // [ratio]
    234293
     294    //
     295    // Concentration around COG (it is calculated here, because the
     296    // distance of the pixel to COG is calculated anyhow)
     297    //
     298    fConcCOG = 0;
     299    for (UInt_t i=0; i<TMath::Min(3U, npix); i++)
     300        fConcCOG += idx[i]<0 ? 0 : evt[idx[i]].GetNumPhotons()*geom.GetPixRatio(idx[i]);
     301    fConcCOG /= hillas.GetSize();                        // [ratio]
     302
     303    // Concentration of signal contained in ellipse
     304    fConcCore /= hillas.GetSize();
     305
     306//    if (fConcCore==0)
     307//        fConcCore = fConcCOG/3;
     308
    235309    SetReadyToSave();
    236 } 
     310}
    237311
    238312// --------------------------------------------------------------------------
    … …  
    249323    *fLog << " - Conc             [1] = " << fConc          << endl;
    250324    *fLog << " - Conc1            [1] = " << fConc1         << endl;
     325    *fLog << " - ConcCOG          [1] = " << fConcCOG       << endl;
     326    *fLog << " - ConcCore         [1] = " << fConcCore      << endl;
    251327    *fLog << " - Num Used Pixels  [#] = " << fNumUsedPixels << endl;
    252328    *fLog << " - Num Core Pixels  [#] = " << fNumCorePixels << endl;
    … …  
    271347    *fLog << " - Conc             [1] = " << fConc          << endl;
    272348    *fLog << " - Conc1            [1] = " << fConc1         << endl;
     349    *fLog << " - ConcCOG          [1] = " << fConcCOG       << endl;
     350    *fLog << " - ConcCore         [1] = " << fConcCore      << endl;
    273351    *fLog << " - Num Used Pixels  [#] = " << fNumUsedPixels << endl;
    274352    *fLog << " - Num Core Pixels  [#] = " << fNumCorePixels << endl;

  • trunk/MagicSoft/Mars/mimage/MNewImagePar.h

    r6855 r7538  
    2323    Float_t fConc;                 // [ratio] concentration ratio: sum of the two highest pixels / fSize
    2424    Float_t fConc1;                // [ratio] concentration ratio: sum of the highest pixel / fSize
     25    Float_t fConcCOG;              // [ratio] concentration of the three pixels next to COG
     26    Float_t fConcCore;             // [ratio] concentration of signals inside or touching the ellipse
    2527
    2628    Float_t fUsedArea;             // Area of pixels which survived the image cleaning
    … …  
    4143    Float_t GetInnerSize()     const { return fInnerSize; }
    4244
    43     Float_t GetConc() const  { return fConc;  }
    44     Float_t GetConc1() const { return fConc1; }
     45    Float_t GetConc() const     { return fConc;     }
     46    Float_t GetConc1() const    { return fConc1;    }
     47    Float_t GetConcCOG() const  { return fConcCOG;  }
     48    Float_t GetConcCore() const { return fConcCore; }
    4549
    4650    Short_t GetNumUsedPixels() const { return fNumUsedPixels; }
    … …  
    5660              const MHillas &hillas, Int_t island=-1);
    5761
    58     ClassDef(MNewImagePar, 4) // Container to hold new image parameters
     62    ClassDef(MNewImagePar, 5) // Container to hold new image parameters
    5963};
    6064
Note: See TracChangeset for help on using the changeset viewer.