Changeset 7178 for trunk/MagicSoft/Mars/mranforest/MRFEnergyEst.cc

Timestamp:

07/08/05 18:14:35 (19 years ago)

Author:

tbretz

Message:

*** empty log message ***

File:

• trunk/MagicSoft/Mars/mranforest/MRFEnergyEst.cc (modified) (12 diffs)

trunk/MagicSoft/Mars/mranforest/MRFEnergyEst.cc

@@ -69 +69 @@
 //
 MRFEnergyEst::MRFEnergyEst(const char *name, const char *title)
-    : fNumTrees(-1), fNumTry(-1), fNdSize(-1), fDebug(kFALSE),
-    fData(0), fEnergyEst(0), fTestMatrix(0), fEstimationMode(kMean)
+    : fDebug(kFALSE), fData(0), fEnergyEst(0),
+    fNumTrees(-1), fNumTry(-1), fNdSize(-1), 
+    fTestMatrix(0), fEstimationMode(kMean)
 {
     fName  = name  ? name  : gsDefName.Data();
     fTitle = title ? title : gsDefTitle.Data();
+}
+
+MRFEnergyEst::~MRFEnergyEst()
+{
+    fEForests.Delete();
 }
 
@@ -163 +169 @@
         MHMatrix matrix1(mat1, "MatrixHadrons");
         MHMatrix matrix0(mat0, "MatrixGammas");
-
-        //matrix1.AddColumns(&usedrules);
-        //matrix0.AddColumns(&usedrules);
 
         // training of RF
@@ -191 +194 @@
             gLog.SetNullOutput(kFALSE);
 
-        const Double_t E = (log10(grid[ie])+log10(grid[ie+1]))/2;
+        // Calculate bin center
+        const Double_t E = (TMath::Log10(grid[ie])+TMath::Log10(grid[ie+1]))/2;
 
         // save whole forest
         MRanForest *forest=(MRanForest*)plist.FindObject("MRanForest");
-        const TString title = Form("%.10f", E);
-        forest->SetTitle(title);
-        forest->Write(title);
+        forest->SetUserVal(E);
+        forest->Write(Form("%.10f", E));
     }
 
@@ -203 +206 @@
     usedrules.Write("rules");
 
-    fileRF.Close();
-
     return kTRUE;
 }
@@ -210 +211 @@
 Int_t MRFEnergyEst::ReadForests(MParList &plist)
 {
-    TFile fileRF(fFileName,"read");
+    TFile fileRF(fFileName, "read");
     if (!fileRF.IsOpen())
     {
@@ -218 +219 @@
 
     fEForests.Delete();
-
-    Int_t i=0;
 
     TIter Next(fileRF.GetListOfKeys());
@@ -230 +229 @@
             continue;
 
-        forest->SetTitle(o->GetTitle());
-        forest->SetBit(kCanDelete);
-
-        fBinning.Set(i+1);
-        fBinning[i++] = atof(o->GetTitle());
+        forest->SetUserVal(atof(o->GetName()));
 
         fEForests.Add(forest);
@@ -283 +278 @@
 //
 //
+#include <TGraph.h>
+#include <TF1.h>
 Int_t MRFEnergyEst::Process()
 {
+    static TF1 f1("f1", "gaus");
+
     TVector event;
     if (fTestMatrix)
@@ -291 +290 @@
         *fData >> event;
 
-    Double_t eest = 0;
-    Double_t hsum = 0;
+    Double_t sume = 0;
+    Double_t sumh = 0;
     Double_t maxh = 0;
     Double_t maxe = 0;
 
-    Int_t i=0;
+    Double_t max  = -1e10;
+    Double_t min  =  1e10;
+
+    //TH1C h("", "", fEForests.GetSize(), 0, 1);
 
     TIter Next(&fEForests);
     MRanForest *rf = 0;
+
+    TGraph g;
     while ((rf=(MRanForest*)Next()))
     {
         const Double_t h = rf->CalcHadroness(event);
-        const Double_t e = fBinning[i++];
-
-        hsum += h;
-        eest += e*h;
+        const Double_t e = rf->GetUserVal();
+        g.SetPoint(g.GetN(), e, h);
+        sume += e*h;
+        sumh += h;
         if (h>maxh)
         {
@@ -312 +316 @@
             maxe = e;
         }
+        if (e>max)
+            max = e;
+        if (e<min)
+            min = e;
     }
 
@@ -317 +325 @@
     {
     case kMean:
-        fEnergyEst->SetVal(pow(10, eest/hsum));
+        fEnergyEst->SetVal(pow(10, sume/sumh));
         break;
     case kMaximum:
         fEnergyEst->SetVal(pow(10, maxe));
         break;
+    case kFit:
+        f1.SetParameter(0, maxh);
+        f1.SetParameter(1, maxe);
+        f1.SetParameter(2, 0.125);
+        g.Fit(&f1, "Q0N");
+        fEnergyEst->SetVal(pow(10, f1.GetParameter(1)));
+        break;
+
     }
     fEnergyEst->SetReadyToSave();
@@ -353 +369 @@
         if (txt==(TString)"maximum")
             fEstimationMode = kMaximum;
+        if (txt==(TString)"fit")
+            fEstimationMode = kFit;
         rc = kTRUE;
     }