Changeset 8679 for trunk/MagicSoft

Timestamp:

08/19/07 22:40:04 (17 years ago)

Author:

tbretz

Message:

*** empty log message ***

Location:

trunk/MagicSoft

Files:

• Mars/Changelog (modified) (1 diff)
• Mars/NEWS (modified) (4 diffs)
• Mars/mhbase/MBinning.cc (modified) (2 diffs)
• Mars/mhbase/MBinning.h (modified) (1 diff)
• Mars/mhbase/MH.cc (modified) (3 diffs)
• Mars/mhbase/MH.h (modified) (1 diff)
• Mars/mjoptim/MJOptimizeDisp.cc (modified) (1 diff)
• Mars/mjoptim/MJOptimizeDisp.h (modified) (1 diff)
• Mars/mjoptim/MJOptimizeEnergy.cc (modified) (4 diffs)
• Mars/mjoptim/MJOptimizeEnergy.h (modified) (1 diff)
• Mars/msignal/MExtractTimeAndChargeSpline.cc (modified) (2 diffs)
• Mars/mtools/MChisqEval.cc (modified) (7 diffs)
• Mars/mtools/MChisqEval.h (modified) (2 diffs)
• include-Classes/MMcFormat/MMcCorsikaRunHeader.cc (modified) (1 diff)

trunk/MagicSoft/Mars/Changelog

@@ -25 +25 @@
        the ListOfCleanups instead of removing the newly created canvas
        from the global list of canvases
+
+   * mimage/MHHillas.cc, mimage/MHHillasExt.cc, mimage/MHHillasSrc.cc,
+     mimage/MHImagePar.cc, mimage/MHNewImagePar.cc:
+     - enables a lot of grids
+
+   * mhbase/MBinning.[h,cc]:
+     - addec constructor with a TArrayD
+
+   * mhbase/MH.[h,cc]:
+     - added two new member functions GetRange and GetRangeUser
+
+   * mjoptim/MJOptimizeDisp.[h,cc]:
+     - added new option UseThetaSq
+
+   * mjoptim/MJOptimizeEnergy.[h,cc]:
+     - added possibility to use also logarithmic fits
+
+   * mmc/MMcCorsikaRunHeader.cc:
+     - improved Print function
+
+   * msignal/MExtractTimeAndChargeSpline.cc:
+     - initialize rise and fall time as suggested by valgrind
+
+   * mtools/MChisqEval.[h,cc]:
+     - implemented the possibility to use weights
+
 
 

trunk/MagicSoft/Mars/NEWS

@@ -132 +132 @@
       + PostCleanType: 3
 
+   - star/ganymed: The old spark-cuts have been replaced by new ones.
+     These new spark cuts have been cross checked with the sequences
+     84720 (before splitter), 101041 (after splitter) and 223539 (new
+     FADCs). In all cases they seperate pretty well, but of course not
+     perfect. Monte Carlos have suggested to tighten the cuts a little
+     bit further at low sizes. This has been done.
+
+   - ganymed: IMPORTANT - The spark cuts have been changed to match
+     all available data more or less well. This might mean that for 
+     YOUR data there are still sparks visible. 
+     Whenever you do an analysis you have to make sure that NO sparks
+     survive your cuts!
+
    - ganymed: The old hadronness and size cuts have gotten new indices
      (10, 11). The old indices 8 and 9 now mean an additional condition
@@ -188 +201 @@
      paremtrization set parameter 8 to 0.
 
-   - star/ganymed: The old spark-cuts have been replaced by new ones.
-     These new spark cuts have been cross checked with the sequences
-     84720 (before splitter), 101041 (after splitter) and 223539 (new
-     FADCs). In all cases they seperate pretty well, but of course not
-     perfect. Monte Carlos have suggested to tighten the cuts a little
-     bit further at low sizes. This has been done.
-
-   - ganymed: IMPORTANT - The spark cuts have been changed to match
-     all available data more or less well. This might mean that for 
-     YOUR data there are still sparks visible. 
-     Whenever you do an analysis you have to make sure that NO sparks
-     survive your cuts!
-
    - ganymed: The distribution of observation time versus zenith
      angle is now shown for off- AND on-data.
@@ -209 +209 @@
      all off-regions from this file. Please do not forget to
      set the scaling in optimwobble.C properly:
-        fit.SetScaleMode(MAlphaFitter::kNone);
         fit.SetScaleUser(1./3);
      it is not yet read in automatically (will follow soon).
@@ -244 +243 @@
      the information is lost.
 
-   - sponde: added a new tab "EventDist" showing the unweighted event
-     distribution of your sample after cuts. This tells you how many
-     events with this energy you had in your MC files. The same information
-     you get from the error bars of the weighted histograms, but this is
-     less intuitive.
+   - sponde: added a new tab "EventDist" showing the unweighted real
+     absolute number of events of your sample after cuts. This tells
+     you how many events with this energy you had in your MC files.
+     The same information you get from the error bars of the weighted
+     histograms, but this is less intuitive.
 
 

trunk/MagicSoft/Mars/mhbase/MBinning.cc

@@ -1 +1 @@
 /* ======================================================================== *\
-! $Name: not supported by cvs2svn $:$Id: MBinning.cc,v 1.17 2006-10-17 17:16:00 tbretz Exp $
+! $Name: not supported by cvs2svn $:$Id: MBinning.cc,v 1.18 2007-08-19 21:35:00 tbretz Exp $
 ! --------------------------------------------------------------------------
 !
@@ -134 +134 @@
 // --------------------------------------------------------------------------
 //
+// Initialize Binning from TArrayD.
+//
+MBinning::MBinning(const TArrayD &axis, const char *name, const char *title)
+{
+    fName  = name  ? name : gsDefName.Data();
+    fTitle = title ? title : gsDefTitle.Data();
+
+    SetEdges(axis);
+}
+
+// --------------------------------------------------------------------------
+//
 // Search in the parameter list for the binning with name "name". If found,
 // set the edges and title accordingly. Default is name of object.

trunk/MagicSoft/Mars/mhbase/MBinning.h

@@ -43 +43 @@
     MBinning(const TH1 &h, const Char_t axis='x', const char *name=0, const char *title=0);
     MBinning(const TAxis &a, const char *name=0, const char *title=0);
+    MBinning(const TArrayD &a, const char *name=0, const char *title=0);
 
     void Copy(TObject &named) const

trunk/MagicSoft/Mars/mhbase/MH.cc

@@ -1 +1 @@
 /* ======================================================================== *\
-! $Name: not supported by cvs2svn $:$Id: MH.cc,v 1.33 2007-08-06 14:44:13 tbretz Exp $
+! $Name: not supported by cvs2svn $:$Id: MH.cc,v 1.34 2007-08-19 21:36:00 tbretz Exp $
 ! --------------------------------------------------------------------------
 !
@@ -20 +20 @@
 !   Author(s): Thomas Bretz  07/2001 <mailto:tbretz@astro.uni-wuerzburg.de>
 !
-!   Copyright: MAGIC Software Development, 2000-2002
+!   Copyright: MAGIC Software Development, 2000-2007
 !
 !
@@ -1409 +1409 @@
 // --------------------------------------------------------------------------
 //
+// Return the first and last bin of the histogram which is not 0
+//
+void MH::GetRange(const TH1 &h, Int_t &lo, Int_t &hi)
+{
+    lo = 0;
+    hi = 1;
+
+    for (int i=2; i<=h.GetNbinsX(); i++)
+    {
+        if (lo==0 && h.GetBinContent(i)>0)
+            lo = i;
+
+        if (h.GetBinContent(i)>0)
+            hi = i;
+    }
+}
+
+// --------------------------------------------------------------------------
+//
+// Return the lower edge of the first and the upper edge of the last bin
+// of the histogram which is not 0
+//
+void MH::GetRangeUser(const TH1 &h, Axis_t &lo, Axis_t &hi)
+{
+    Int_t f, l;
+    GetRange(h, f, l);
+
+    lo = h.GetBinLowEdge(f);
+    hi = h.GetBinLowEdge(l+1);
+}
+
+// --------------------------------------------------------------------------
+//
 // See MTask::PrintSkipped
 //

trunk/MagicSoft/Mars/mhbase/MH.h

@@ -119 +119 @@
     static void ProjectionY(TH1D &dest, const TH2 &src, Int_t firstxbin=-1, Int_t lastxbin=9999);
 
+    static void GetRangeUser(const TH1 &hist, Axis_t &lo, Axis_t &hi);
+    static void GetRange(const TH1 &hist, Int_t &lo, Int_t &hi);
+
     static TObject *FindObjectInPad(const char *name, TVirtualPad *pad=NULL);
 

trunk/MagicSoft/Mars/mjoptim/MJOptimizeDisp.cc

@@ -178 +178 @@
 
     MChisqEval eval;
-    eval.SetY1("sqrt(ThetaSquared.fVal)");
+    eval.SetY1(fUseThetaSq?"ThetaSquared.fVal":"sqrt(ThetaSquared.fVal)");
+    if (weights)
+        eval.SetNameWeight();
 
     MMatrixLoop loop(&m);

trunk/MagicSoft/Mars/mjoptim/MJOptimizeDisp.h

@@ -10 +10 @@
 class MJOptimizeDisp : public MJOptimize
 {
+private:
+    Bool_t fUseThetaSq;
+
 public:
-    MJOptimizeDisp() : MJOptimize() { }
+    MJOptimizeDisp() : MJOptimize(), fUseThetaSq(kFALSE) { }
+
+    void EnableThetaSq(Bool_t b=kTRUE) { fUseThetaSq=b; }
 
     // Special optimizing routines

trunk/MagicSoft/Mars/mjoptim/MJOptimizeEnergy.cc

@@ -66 +66 @@
 // histograms
 #include "../mhflux/MHEnergyEst.h"
+#include "../mtools/MChisqEval.h"
 
 // tasks
@@ -108 +109 @@
     MHMatrix m("M");
     AddRulesToMatrix(m);
+    const Int_t map = m.AddColumn("MMcEvt.fEnergy");
     parlist.AddToList(&m);
 
@@ -136 +138 @@
         fill.SetWeight();
 
+    MChisqEval eval;
+    eval.SetY1(fOptimLog?Form("log10(MEnergyEst.fVal/M[%d])", map):Form("MEnergyEst.fVal-M[%d]", map));
+    if (weights)
+        eval.SetNameWeight();
+
     MMatrixLoop loop(&m);
 
@@ -143 +150 @@
         tasklist.AddToList(weights);
     tasklist.AddToList(&fill);
+    tasklist.AddToList(&eval);
 
     // Optimize with the tasklist in this parameterlist

trunk/MagicSoft/Mars/mjoptim/MJOptimizeEnergy.h

@@ -10 +10 @@
 class MJOptimizeEnergy : public MJOptimize
 {
+private:
+    Bool_t fOptimLog;
+
 public:
-    MJOptimizeEnergy() : MJOptimize() { }
+    MJOptimizeEnergy() : MJOptimize(), fOptimLog(kFALSE) { }
+
+    void EnableOptimLog(Bool_t b=kTRUE) { fOptimLog=b; }
 
     // Special optimizing routines

trunk/MagicSoft/Mars/msignal/MExtractTimeAndChargeSpline.cc

@@ -1 +1 @@
 /* ======================================================================== *\
-! $Name: not supported by cvs2svn $:$Id: MExtractTimeAndChargeSpline.cc,v 1.68 2007-05-09 12:15:53 tbretz Exp $
+! $Name: not supported by cvs2svn $:$Id: MExtractTimeAndChargeSpline.cc,v 1.69 2007-08-19 21:40:03 tbretz Exp $
 ! --------------------------------------------------------------------------
 !
@@ -179 +179 @@
 //
 MExtractTimeAndChargeSpline::MExtractTimeAndChargeSpline(const char *name, const char *title) 
-    : fHeightTm(0.5), fExtractionType(MExtralgoSpline::kIntegralRel)
+    : fRiseTimeHiGain(0), fFallTimeHiGain(0), fHeightTm(0.5), fExtractionType(MExtralgoSpline::kIntegralRel)
 {
 

trunk/MagicSoft/Mars/mtools/MChisqEval.cc

@@ -49 +49 @@
 const TString MChisqEval::gsDefTitle = "Evaluate a chisq";
 
-MChisqEval::MChisqEval(const char *name, const char *title) : fData0(NULL), fData1(NULL), fNameResult("MinimizationValue")
+MChisqEval::MChisqEval(const char *name, const char *title)
+    : fData0(NULL), fData1(NULL), fWeight(NULL), fNameResult("MinimizationValue")
 {
     fName  = name  ? name  : gsDefName.Data();
@@ -55 +56 @@
 }
 
-MChisqEval::MChisqEval(MData *y1, const char *name, const char *title) : fData0(NULL), fData1(NULL), fNameResult("MinimizationValue")
+MChisqEval::MChisqEval(MData *y1, const char *name, const char *title)
+: fData0(NULL), fData1(NULL), fWeight(NULL), fNameResult("MinimizationValue")
 {
     fName  = name  ? name  : gsDefName.Data();
@@ -62 +64 @@
 }
 
-MChisqEval::MChisqEval(MData *y1, MData *y2, const char *name, const char *title) : fData0(NULL), fData1(NULL), fNameResult("MinimizationValue")
+MChisqEval::MChisqEval(MData *y1, MData *y2, const char *name, const char *title)
+: fData0(NULL), fData1(NULL), fWeight(NULL), fNameResult("MinimizationValue")
 {
     fName  = name  ? name  : gsDefName.Data();
@@ -112 +115 @@
 {
     fChisq = 0;
+    fSumW  = 0;
 
     if (!fData0)
@@ -122 +126 @@
         if (!fData1->PreProcess(plist))
             return kFALSE;
+
+    if (!fNameWeight.IsNull())
+    {
+        fWeight = (MParameterD*)plist->FindObject(fNameWeight, "MParameterD");
+        if (!fWeight)
+            return kFALSE;
+    }
 
     fResult = (MParameterD*)plist->FindCreateObj("MParameterD", fNameResult);
@@ -138 +149 @@
     const Double_t err = fData1 ? y1*y2 : 1;
 
-    fChisq += dy*dy/err;
+    const Double_t w   = fWeight ? fWeight->GetVal() : 1;
+
+    fChisq += w*dy*dy/err;
+    fSumW  += w;
+
     return kTRUE;
 }
@@ -145 +160 @@
 {
     if (GetNumExecutions()>0)
-        fChisq /= GetNumExecutions();
+        fChisq /= fSumW;
 
     fResult->SetVal(fChisq);

trunk/MagicSoft/Mars/mtools/MChisqEval.h

@@ -15 +15 @@
     static const TString gsDefTitle;
 
-    Double_t     fChisq;  //! Evaluated chi square
-    MParameterD *fResult; //! Storage for result
-
     MData   *fData0; // Data Member one (monte carlo data or chisq function)
     MData   *fData1; // Data Member two (measured data)
 
+    MParameterD *fWeight;  //! Storage for weight
+    MParameterD *fResult;  //! Storage for result
+
     TString fNameResult;
+    TString fNameWeight;
+
+    Double_t     fChisq;   //! Evaluated chi square
+    Double_t     fSumW;    //! Sum of weights
 
     void StreamPrimitive(ostream &out) const;
@@ -43 +47 @@
 
     void SetOwner(Bool_t o=kTRUE) { o ? SetBit(kIsOwner) : ResetBit(kIsOwner); }
+    void SetNameWeight(TString w="MWeight") { fNameWeight=w; }
 
     Double_t GetChisq() const { return fChisq; }

trunk/MagicSoft/include-Classes/MMcFormat/MMcCorsikaRunHeader.cc

@@ -206 +206 @@
     //   Print Information about the Geometry of the camera
     //
-  *fLog << all << GetTitle() <<". Number of Telescopes: " << fNumCT << endl;
-
-  fTelescopes.Print();
+    *fLog << all << GetTitle() <<":" << endl;
+    *fLog << " Spectral Slope: " << fSlopeSpec << " from " << fELowLim << "GeV to " << fEUppLim << "GeV" << endl;
+    *fLog << " Number of Telescopes: " << fNumCT << endl;
+    fTelescopes.Print();
 }