Changeset 7710 for trunk/MagicSoft

Timestamp:

05/18/06 14:01:43 (18 years ago)

Author:

tbretz

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars

Files:

• Changelog (modified) (1 diff)
• mhflux/MAlphaFitter.h (modified) (1 diff)
• mhflux/MHThetaSq.cc (modified) (1 diff)
• mhflux/MHThetaSq.h (modified) (1 diff)
• mjtrain/MJTrainRanForest.cc (modified) (1 diff)
• mranforest/MRanForestCalc.cc (modified) (4 diffs)
• mranforest/MRanTree.cc (modified) (11 diffs)

trunk/MagicSoft/Mars/Changelog

@@ -18 +18 @@
 
                                                  -*-*- END OF LINE -*-*-
+
+ 2006/05/18 Thomas Bretz
+
+   * mhflux/MAlphaFitter.h:
+     - added Getter for ScaleMode
+
+   * mhflux/MHThetaSq.[h,cc]:
+     - removed obsolete data-member fThetaSq
+
+   * mjtrain/MJTrainRanForest.cc, mranforest/MRanForestCalc.cc:
+     - fixed a typo in a comment
+
+   * mranforest/MRanTree.cc:
+     - added a lot of comments in the code
+
+
 
  2006/05/18 Daniela Dorner

trunk/MagicSoft/Mars/mhflux/MAlphaFitter.h

@@ -169 +169 @@
     Double_t GetMinimizationValue() const;
 
+    ScaleMode_t GetScaleMode() const       { return fScaleMode; }
+
     Double_t GetGausSigma() const          { return fCoefficients[2]; }
     Double_t GetGausMu() const             { return fCoefficients[1]; }

trunk/MagicSoft/Mars/mhflux/MHThetaSq.cc

@@ -62 +62 @@
 //
 MHThetaSq::MHThetaSq(const char *name, const char *title)
-    : MHAlpha(name, title), fThetaSq(0), fNumBinsSignal(3), fNumBinsTotal(75)
+    : MHAlpha(name, title), fNumBinsSignal(3), fNumBinsTotal(75)
 {
     //

trunk/MagicSoft/Mars/mhflux/MHThetaSq.h

@@ -11 +11 @@
 {
 private:
-    MParameterD  *fThetaSq; //!
-
     UInt_t fNumBinsSignal;
     UInt_t fNumBinsTotal;

trunk/MagicSoft/Mars/mjtrain/MJTrainRanForest.cc

@@ -90 +90 @@
 //   Int_t idx = AddParameter("log10(MHillas.fSize)");
 //
-// The indices area starting with 0 always.
+// The indices are starting with 0 always.
 //
 Int_t MJTrainRanForest::AddParameter(const char *rule)

trunk/MagicSoft/Mars/mranforest/MRanForestCalc.cc

@@ -1 @@
-/* ======================================================================== *\
+/* ================================q======================================== *\
 !
 ! *
@@ -32 +32 @@
 #include "MRanForestCalc.h"
 
+#include <TF1.h>
+#include <TGraph.h>
 #include <TVector.h>
 
@@ -279 +281 @@
     }
 
-    // Maybe fEForests[0].fRules yould be used instead?
-
+    // Maybe fEForests[0].fRules could be used instead?
     if (fData->Read("rules")<=0)
     {
@@ -322 +323 @@
 }
 
-#include <TGraph.h>
-#include <TF1.h>
 Int_t MRanForestCalc::Process()
 {

trunk/MagicSoft/Mars/mranforest/MRanTree.cc

@@ -41 +41 @@
 #include "MArrayF.h"
 
+#include "MMath.h"
+
 #include "MLog.h"
 #include "MLogManip.h"
@@ -176 +178 @@
     Double_t pdo=0;
 
-    for (Int_t j=0; j<nclass; j++)
+    // tclasspop: sum of weights for events in class
+    for (Int_t j=0; j<nclass; j++) // loop over number of classes to classifiy
     {
         pno+=tclasspop[j]*tclasspop[j];
@@ -182 +185 @@
     }
 
-    const Double_t crit0=pno/pdo;
+    const Double_t crit0=pno/pdo;  // weighted mean of weights
 
     // start main loop through variables to find best split,
@@ -190 +193 @@
 
     // random split selection, number of trials = fNumTry
-    for (Int_t mt=0; mt<fNumTry; mt++)
+    for (Int_t mt=0; mt<fNumTry; mt++) // we could try ALL variables???
     {
         const Int_t mvar= gRandom->Integer(mdim);
@@ -212 +215 @@
 
             // do classification, Gini index as split rule
-            rln+=u*( 2*wl[k]+u);
-            rrn+=u*(-2*wr[k]+u);
-
-            rld+=u;
-            rrd-=u;
-
-            wl[k]+=u;
-            wr[k]-=u;
-
+            rln   +=u*(2*wl[k]+u);  // += u*(wl[k]{i-1} + wl[k]{i-1}+u{i})
+            rld   +=u;   // sum of weights left  from cut total
+            wl[k] +=u;   // sum of weights left  from cut for class k
+
+            rrn   -=u*(2*wr[k]-u);  // -= u*(wr[k]{i-1} + wr[k]{i-1}-u{i})
+            //  rr0=0; rr0+=u*2*tclasspop[k]
+            //  rrn = pno - rr0 + rln
+            rrd   -=u;   // sum of weights right from cut total
+            wr[k] -=u;   // sum of weights right from cut for class k
+
+            // REPLACE BY?
+            // rr0   = 0
+            // rr0  += u*2*tclasspop[k]
+            // rrn   = pno - rr0 + rln
+            // rrd   = pdo - rld
+            // wr[k] = tclasspop[k] - wl[k]
+
+            // crit = (rln*(pdo - rld + 1) + pno - rr0) / rld*(pdo - rld)
+
+            /*
+             if (k==background)
+                continue;
+             crit = TMath::Max(MMath::SignificanceLiMa(rld, rld-wl[k]),
+                               MMath::SignificanceLiMa(rrd, rrd-wr[k]))
+             */
+
+            // This condition is in fact a == (> cannot happen at all)
+            // This is because we cannot set the cut between two identical values
+            //if (datarang[mn+datasort[mn+nsp]]>=datarang[mn+datasort[mn+nsp+1]])
             if (datarang[mn+nc]>=datarang[mn+datasort[mn+nsp+1]])
                 continue;
 
-            if (TMath::Min(rrd,rld)<=1.0e-5)
+            // If crit starts to become pretty large do WHAT???
+            if (TMath::Min(rrd,rld)<=1.0e-5) // FIXME: CHECKIT FOR WEIGHTS!
                 continue;
 
             const Double_t crit=(rln/rld)+(rrn/rrd);
 
+            // Search for the highest value of crit
             if (crit<=critvar) continue;
 
+            // store the highest crit value and the corresponding event to cut at
             nbestvar=nsp;
             critvar=crit;
@@ -237 +263 @@
         if (critvar<=critmax) continue;
 
-        msplit=mvar;
-        nbest=nbestvar;
+        msplit=mvar;      // Variable in which to split
+        nbest=nbestvar;   // event at which the best split was found
         critmax=critvar;
     }
 
+    // crit0 = MMath::SignificanceLiMa(pdo, pdo-tclasspop[0])
+    // mean increase of sensitivity
+    // decsplit = sqrt(critmax/crit0)
     decsplit=critmax-crit0;
 
@@ -299 +328 @@
         const Int_t mn  = mvar*numdata;
 
-        Double_t rrn=0, rrd=0, rln=0, rld=0;
-
         Double_t esumr =mean;
         Double_t e2sumr=square;
@@ -316 +343 @@
             const Float_t &u=winbag[nc];
 
-            e2sumr-=u*f*f;
-            esumr -=u*f;
-            wr    -=u;
-
-            //-------------------------------------------
-            // resolution
-            //rrn=(wr*e2sumr-esumr*esumr)*wr;
-            //rrd=(wr-1.)*esumr*esumr;
-
-            // resolution times n
-            //rrn=(wr*e2sumr-esumr*esumr)*wr;
-            //rrd=esumr*esumr;
-
-            // sigma
-            //rrn=(e2sumr-esumr*esumr/wr);
-            //rrd=(wr-1.);
-
-            // sigma times n
-            rrn=(e2sumr-esumr*esumr/wr);
-            rrd=1.;
-
-            // 1./(n*variance)
-            //rrn=1.;
-            //rrd=(e2sumr-esumr*esumr/wr);
-            //-------------------------------------------
-
             e2suml+=u*f*f;
             esuml +=u*f;
@@ -348 +349 @@
             //-------------------------------------------
             // resolution
-            //rln=(wl*e2suml-esuml*esuml)*wl;
-            //rld=(wl-1.)*esuml*esuml;
+            //const Double_t rln=(wl*e2suml-esuml*esuml)*wl;
+            //const Double_t rld=(wl-1.)*esuml*esuml;
 
             // resolution times n
-            //rln=(wl*e2suml-esuml*esuml)*wl;
-            //rld=esuml*esuml;
+            //const Double_t rln=(wl*e2suml-esuml*esuml)*wl;
+            //const Double_t rld=esuml*esuml;
 
             // sigma
-            //rln=(e2suml-esuml*esuml/wl);
-            //rld=(wl-1.);
+            //const Double_t rln=(e2suml-esuml*esuml/wl);
+            //const Double_t rld=(wl-1.);
 
             // sigma times n
-            rln=(e2suml-esuml*esuml/wl);
-            rld=1.;
+            Double_t rln=(e2suml-esuml*esuml/wl);
+            Double_t rld=1.;
 
             // 1./(n*variance)
-            //rln=1.;
-            //rld=(e2suml-esuml*esuml/wl);
+            //const Double_t rln=1.;
+            //const Double_t rld=(e2suml-esuml*esuml/wl);
+            //-------------------------------------------
+
+            // REPLACE BY??? 
+            e2sumr-=u*f*f;   // e2sumr = square       - e2suml
+            esumr -=u*f;     // esumr  = mean         - esuml
+            wr    -=u;       // wr     = tclasspop[0] - wl
+
+            //-------------------------------------------
+            // resolution
+            //const Double_t rrn=(wr*e2sumr-esumr*esumr)*wr;
+            //const Double_t rrd=(wr-1.)*esumr*esumr;
+
+            // resolution times n
+            //const Double_t rrn=(wr*e2sumr-esumr*esumr)*wr;
+            //const Double_t rrd=esumr*esumr;
+
+            // sigma
+            //const Double_t rrn=(e2sumr-esumr*esumr/wr);
+            //const Double_t rrd=(wr-1.);
+
+            // sigma times n
+            const Double_t rrn=(e2sumr-esumr*esumr/wr);
+            const Double_t rrd=1.;
+
+            // 1./(n*variance)
+            //const Double_t rrn=1.;
+            //const Double_t rrd=(e2sumr-esumr*esumr/wr);
             //-------------------------------------------
 
@@ -538 +566 @@
               const int j=idclass[nc];
                    
+              // statistics left from cut
               mean[ncur+1]+=hadtrue[nc]*winbag[nc];
               square[ncur+1]+=hadtrue[nc]*hadtrue[nc]*winbag[nc];
 
+              // sum of weights left from cut
               classpop[j*nrnodes+ncur+1]+=winbag[nc];
           }
@@ -549 +579 @@
               const int j=idclass[nc];
 
+              // statistics right from cut
               mean[ncur+2]  +=hadtrue[nc]*winbag[nc];
               square[ncur+2]+=hadtrue[nc]*hadtrue[nc]*winbag[nc];
 
+              // sum of weights right from cut
               classpop[j*nrnodes+ncur+2]+=winbag[nc];
           }