Index: trunk/MagicSoft/Mars/mimage/MHNewImagePar.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHNewImagePar.cc (revision 7535) +++ trunk/MagicSoft/Mars/mimage/MHNewImagePar.cc (revision 7538) @@ -19,5 +19,5 @@ ! Author(s): Thomas Bretz, 04/2003 ! -! Copyright: MAGIC Software Development, 2000-2004 +! Copyright: MAGIC Software Development, 2000-2006 ! ! @@ -27,4 +27,10 @@ // // MHNewImagePar +// +// +// Version 2: +// ---------- +// + TH1F fHistConcCOG; // [ratio] concentration around the center of gravity (all pixels) +// + TH1F fHistConcCore; // [ratio] concentration of signals inside ellipse (used pixels) // //////////////////////////////////////////////////////////////////////////// @@ -116,18 +122,32 @@ fHistConc.SetDirectory(NULL); fHistConc1.SetDirectory(NULL); + fHistConcCOG.SetDirectory(NULL); + fHistConcCore.SetDirectory(NULL); fHistConc.SetName("Conc2"); fHistConc1.SetName("Conc1"); + fHistConcCOG.SetName("ConcCOG"); + fHistConcCore.SetName("ConcCore"); fHistConc.SetTitle("Ratio: Conc"); fHistConc1.SetTitle("Ratio: Conc1"); + fHistConcCOG.SetTitle("Ratio: ConcCOG"); + fHistConcCore.SetTitle("Ratio: ConcCore"); fHistConc.SetXTitle("Ratio"); fHistConc1.SetXTitle("Ratio"); + fHistConcCOG.SetXTitle("Ratio"); + fHistConcCore.SetXTitle("Ratio"); fHistConc.SetYTitle("Counts"); fHistConc1.SetYTitle("Counts"); + fHistConcCOG.SetYTitle("Counts"); + fHistConcCore.SetYTitle("Counts"); + fHistConc.UseCurrentStyle(); fHistConc1.UseCurrentStyle(); - fHistConc.UseCurrentStyle(); + fHistConcCOG.UseCurrentStyle(); + fHistConcCore.UseCurrentStyle(); fHistConc.SetFillStyle(4000); fHistConc1.SetFillStyle(4000); + fHistConcCOG.SetFillStyle(4000); + fHistConcCore.SetFillStyle(4000); fHistConc1.SetLineColor(kBlue); - + fHistConcCOG.SetLineColor(kBlue); MBinning bins; @@ -138,4 +158,6 @@ bins.Apply(fHistConc); bins.Apply(fHistConc1); + bins.Apply(fHistConcCOG); + bins.Apply(fHistConcCore); bins.SetEdges(75, 0.5, 150.5); @@ -189,6 +211,8 @@ //ApplyBinning(*plist, "Area", &fHistCoreArea); - ApplyBinning(*plist, "Conc", &fHistConc); - ApplyBinning(*plist, "Conc1", &fHistConc1); + ApplyBinning(*plist, "Conc", &fHistConc); + ApplyBinning(*plist, "Conc1", &fHistConc1); + ApplyBinning(*plist, "ConcCOG", &fHistConcCOG); + ApplyBinning(*plist, "ConcCore", &fHistConcCore); return kTRUE; @@ -223,4 +247,6 @@ fHistConc.Fill(h.GetConc(), w); fHistConc1.Fill(h.GetConc1(), w); + fHistConcCOG.Fill(h.GetConcCOG(), w); + fHistConcCore.Fill(h.GetConcCore(), w); return kTRUE; @@ -284,6 +310,6 @@ void MHNewImagePar::Paint(Option_t *o) { - if (fHistLeakage1.GetMaximum()>0 && gPad->GetPad(1)) - gPad->GetPad(1)->SetLogy(); + if (fHistLeakage1.GetMaximum()>0 && gPad->GetPad(1) && gPad->GetPad(1)->GetPad(1)) + gPad->GetPad(1)->GetPad(1)->SetLogy(); } @@ -306,5 +332,5 @@ if (!same) - pad->Divide(2,2); + pad->Divide(2, 1); else { @@ -315,4 +341,6 @@ fHistUsedArea.SetName("UsedAreaSame"); fHistCoreArea.SetName("CoreAreaSame"); + fHistConcCOG.SetName("ConcCOGSame"); + fHistConcCore.SetName("ConcCoreSame"); fHistConc1.SetName("Conc1Same"); fHistConc.SetName("Conc2Same"); @@ -324,4 +352,6 @@ fHistUsedArea.SetDirectory(0); fHistCoreArea.SetDirectory(0); + fHistConcCOG.SetDirectory(0); + fHistConcCore.SetDirectory(0); fHistConc1.SetDirectory(0); fHistConc.SetDirectory(0); @@ -332,4 +362,6 @@ fHistUsedPix.SetLineColor(kCyan); fHistConc1.SetLineColor(kMagenta); + fHistConcCOG.SetLineColor(kCyan); + fHistConcCore.SetLineColor(kCyan); fHistConc.SetLineColor(kCyan); fHistCoreArea.SetLineColor(kMagenta); @@ -338,5 +370,11 @@ pad->cd(1); + TVirtualPad *pad1=gPad; + pad1->SetBorderMode(0); + pad1->Divide(1,3, 0.001, 0.001); + + pad1->cd(1); gPad->SetBorderMode(0); + gPad->SetGridx(); TAxis &x = *fHistLeakage1.GetXaxis(); x.SetRangeUser(0.0, x.GetXmax()); @@ -348,21 +386,36 @@ fHistLeakage2.SetMaximum(0.1); // dummy value to allow log-scale - pad->cd(2); + pad1->cd(2); gPad->SetBorderMode(0); + gPad->SetGridx(); RemoveFromPad("UsedPixSame"); RemoveFromPad("CorePixSame"); MH::DrawSame(fHistCorePix, fHistUsedPix, "Number of core/used Pixels", same); - pad->cd(3); + pad1->cd(3); gPad->SetBorderMode(0); + gPad->SetGridx(); + RemoveFromPad("CoreAreaSame"); + RemoveFromPad("UsedAreaSame"); + MH::DrawSame(fHistCoreArea, fHistUsedArea, "Area of core/used Pixels", same); + + pad->cd(2); + TVirtualPad *pad2=gPad; + pad2->SetBorderMode(0); + pad2->Divide(1, 2, 0.001, 0.001); + + pad2->cd(1); + gPad->SetBorderMode(0); + gPad->SetGridx(); RemoveFromPad("Conc1Same"); RemoveFromPad("Conc2Same"); MH::DrawSame(fHistConc1, fHistConc, "Concentrations", same); - pad->cd(4); + pad2->cd(2); gPad->SetBorderMode(0); - RemoveFromPad("CoreAreaSame"); - RemoveFromPad("UsedAreaSame"); - MH::DrawSame(fHistCoreArea, fHistUsedArea, "Area of core/used Pixels", same); + gPad->SetGridx(); + RemoveFromPad("ConcCOGSame"); + RemoveFromPad("ConcCoreSame"); + MH::DrawSame(fHistConcCore, fHistConcCOG, "Concentrations", same); } @@ -373,4 +426,8 @@ if (name.Contains("Leakage2", TString::kIgnoreCase)) return const_cast(&fHistLeakage2); + if (name.Contains("ConcCOG", TString::kIgnoreCase)) + return const_cast(&fHistConcCOG); + if (name.Contains("ConcCore", TString::kIgnoreCase)) + return const_cast(&fHistConcCore); if (name.Contains("Conc1", TString::kIgnoreCase)) // must be first! return const_cast(&fHistConc1); Index: trunk/MagicSoft/Mars/mimage/MHNewImagePar.h =================================================================== --- trunk/MagicSoft/Mars/mimage/MHNewImagePar.h (revision 7535) +++ trunk/MagicSoft/Mars/mimage/MHNewImagePar.h (revision 7538) @@ -25,4 +25,6 @@ TH1F fHistConc; // [ratio] concentration ratio: sum of the two highest pixels / fSize TH1F fHistConc1; // [ratio] concentration ratio: sum of the highest pixel / fSize + TH1F fHistConcCOG; // [ratio] concentration of the three pixels next to COG + TH1F fHistConcCore; // [ratio] concentration of signals inside or touching the ellipse Float_t fMm2Deg; @@ -53,4 +55,6 @@ TH1F &GetHistConc() { return fHistConc; } TH1F &GetHistConc1() { return fHistConc1; } + TH1F &GetHistConcCOG() { return fHistConcCOG; } + TH1F &GetHistConcCore() { return fHistConcCore; } void SetMmScale(Bool_t mmscale=kTRUE); @@ -60,5 +64,5 @@ void Paint(Option_t *opt=""); - ClassDef(MHNewImagePar, 1) // Histograms of new image parameters + ClassDef(MHNewImagePar, 2) // Histograms of new image parameters }; Index: trunk/MagicSoft/Mars/mimage/MHillasExt.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MHillasExt.cc (revision 7535) +++ trunk/MagicSoft/Mars/mimage/MHillasExt.cc (revision 7538) @@ -161,5 +161,5 @@ const Double_t dist = dx*dx+dy*dy; - if (TMath::Abs(dist)>TMath::Abs(maxdist)) + if (dist>TMath::Abs(maxdist)) maxdist = dzx<0 ? -dist : dist; // [mm^2] @@ -177,5 +177,4 @@ maxpix = nphot; // [1] maxpixid = i;//pixid; - continue; // [1] } } @@ -183,4 +182,7 @@ const MGeomPix &maxp = geom[maxpixid]; + // + // Asymmetry + // fAsym = (hil.GetMeanX()-maxp.GetX())*hil.GetCosDelta() + (hil.GetMeanY()-maxp.GetY())*hil.GetSinDelta(); // [mm] @@ -195,4 +197,7 @@ fM3Trans = m3y<0 ? -pow(-m3y, 1./3) : pow(m3y, 1./3); // [mm] + // + // Distance between max signal and COG + // const Double_t md = TMath::Sqrt(TMath::Abs(maxdist)); fMaxDist = maxdist<0 ? -md : md; // [mm] Index: trunk/MagicSoft/Mars/mimage/MNewImagePar.cc =================================================================== --- trunk/MagicSoft/Mars/mimage/MNewImagePar.cc (revision 7535) +++ trunk/MagicSoft/Mars/mimage/MNewImagePar.cc (revision 7538) @@ -38,4 +38,6 @@ // Float_t fConc; // [ratio] concentration ratio: sum of the two highest pixels / fSize // Float_t fConc1; // [ratio] concentration ratio: sum of the highest pixel / fSize +// Float_t fConcCOG; // [ratio] concentration of the three pixels next to COG +// Float_t fConcCore; // [ratio] concentration of signals inside or touching the ellipse // // Float_t fUsedArea; // Area of pixels which survived the image cleaning @@ -65,4 +67,9 @@ // + removed fNumSaturatedPixels // +// Version 5: +// ---------- +// - added fConcCOG +// - added fConcCore +// // ///////////////////////////////////////////////////////////////////////////// @@ -105,8 +112,10 @@ fInnerLeakage1 = -1; fInnerLeakage2 = -1; - fInnerSize = -1; - - fConc = -1; - fConc1 = -1; + fInnerSize = -1; + + fConc = -1; + fConc1 = -1; + fConcCOG = -1; + fConcCore = -1; fNumUsedPixels = -1; @@ -131,4 +140,5 @@ fInnerSize = 0; + fConcCore = 0; Double_t edgepix1 = 0; @@ -140,4 +150,11 @@ Float_t maxpix1 = 0; // [#phot] Float_t maxpix2 = 0; // [#phot] + + const Double_t d = geom.GetMaxRadius()*geom.GetMaxRadius(); + Double_t dist[3] = { d, d, d }; + Int_t idx[3] = { -1, -1, -1}; + + const Double_t rl = 1./(hillas.GetLength()*hillas.GetLength()); + const Double_t rw = 1./(hillas.GetWidth() *hillas.GetWidth()); const Bool_t ismagiclike = @@ -149,4 +166,39 @@ for (UInt_t i=0; i=0 && pix.GetIdxIsland()!=island) continue; - - // Get geometry of pixel - //const Int_t pixid = pix->GetPixId(); - const MGeomPix &gpix = geom[i/*pixid*/]; // count used and core pixels @@ -172,5 +220,17 @@ fUsedArea += gpix.GetA(); + // signal in pixel Double_t nphot = pix.GetNumPhotons(); + + // + // Calculate signal contained inside ellipse + // + const Double_t dzx = hillas.GetCosDelta()*dx + hillas.GetSinDelta()*dy; // [mm] + const Double_t dzy = -hillas.GetSinDelta()*dx + hillas.GetCosDelta()*dy; // [mm] + const Double_t dz = gpix.GetD()*gpix.GetD()/4; + const Double_t tana = dzy*dzy/(dzx*dzx); + const Double_t distr = (1+tana)/(rl + tana*rw); + if (distr>dist0-dz || dzx==0) + fConcCore += nphot; // @@ -213,11 +273,10 @@ maxpix2 = maxpix1; maxpix1 = nphot; // [1] - continue; // [1] } - - if (nphot>maxpix2) - maxpix2 = nphot; // [1] + else + if (nphot>maxpix2) + maxpix2 = nphot; // [1] } - + fInnerLeakage1 = edgepixin1 / fInnerSize; fInnerLeakage2 = edgepixin2 / fInnerSize; @@ -226,13 +285,28 @@ fLeakage2 = edgepix2 / hillas.GetSize(); - // FIXME?: in case the pixel with highest signal density is an outer pixel, - // the value of fConc (ratio of signal in two highest pixels to SIZE) should - // rather be 2*fConc1, under the simplest assumption that the light density - // inside the outer (large) pixel is uniform. + // FIXME?: in case the pixel with highest signal density is an outer + // pixel, the value of fConc (ratio of signal in two highest pixels + // to SIZE) should rather be 2*fConc1, under the simplest assumption + // that the light density inside the outer (large) pixel is uniform. fConc = (maxpix1+maxpix2)/hillas.GetSize(); // [ratio] fConc1 = maxpix1/hillas.GetSize(); // [ratio] + // + // Concentration around COG (it is calculated here, because the + // distance of the pixel to COG is calculated anyhow) + // + fConcCOG = 0; + for (UInt_t i=0; i