Changeset 1336

Timestamp:

06/03/02 09:57:00 (23 years ago)

Author:

tbretz

Message:

*** empty log message ***

Location:

trunk/MagicSoft/Mars/mhist

Files:

• HistLinkDef.h (modified) (1 diff)
• MFillH.cc (modified) (1 diff)
• MH3.cc (modified) (4 diffs)
• MHCompProb.cc (added)
• MHCompProb.h (added)
• MHHadroness.cc (added)
• MHHadroness.h (added)
• MHMatrix.cc (modified) (7 diffs)
• MHMatrix.h (modified) (3 diffs)
• Makefile (modified) (2 diffs)

trunk/MagicSoft/Mars/mhist/HistLinkDef.h

@@ -41 +41 @@
 #pragma link C++ class MHMcEnergyMigration+;
 
+#pragma link C++ class MHCompProb+;
+#pragma link C++ class MHHadroness+;
 
 #endif

trunk/MagicSoft/Mars/mhist/MFillH.cc

@@ -243 +243 @@
 
         if (!obj)
+        {
+            if (cls==name)
+                *fLog << inf << "Object '" << fHName << "' not found in parlist... trying to create one." << endl;
             obj = pList->FindCreateObj(cls, name);
+        }
 
         if (!obj)

trunk/MagicSoft/Mars/mhist/MH3.cc

@@ -319 +319 @@
         p->Draw("same");
         p->SetBit(kCanDelete);
+        p->SetDirectory(NULL);
     }
     if (str.Contains("PROFY", TString::kIgnoreCase) && fDimension==2)
@@ -325 +326 @@
         p->Draw("same");
         p->SetBit(kCanDelete);
+        p->SetDirectory(NULL);
     }
 
@@ -352 +354 @@
         p->Draw("same");
         p->SetBit(kCanDelete);
+        p->SetDirectory(NULL);
     }
     if (str.Contains("PROFY", TString::kIgnoreCase) && fDimension==2)
@@ -358 +361 @@
         p->Draw("same");
         p->SetBit(kCanDelete);
+        p->SetDirectory(NULL);
     }
 

trunk/MagicSoft/Mars/mhist/MHMatrix.cc

@@ -43 +43 @@
 //
 /////////////////////////////////////////////////////////////////////////////
-
 #include "MHMatrix.h"
 
 #include <TList.h>
+#include <TArrayD.h>
+#include <TArrayI.h>
 
 #include "MLog.h"
@@ -61 +62 @@
 // Default Constructor
 //
-MHMatrix::MHMatrix()
+MHMatrix::MHMatrix(const char *name, const char *title)
     : fNumRow(0)
 {
+    fName  = name  ? name  : "MHMatrix";
+    fTitle = title ? title : "Multidimensional Matrix";
+
     fData = new TList;
     fData->SetOwner();
+
+    fRules = new TList;
+    fRules->SetOwner();
 }
 
@@ -75 +82 @@
 {
     delete fData;
+    delete fRules;
 }
 
@@ -91 +99 @@
     }
 
-    MDataChain &chain = *new MDataChain(name);
+    MDataChain &chain = *new MDataChain(rule);
     if (!chain.IsValid())
     {
@@ -100 +108 @@
 
     fData->Add(&chain);
+
+    TNamed *name = new TNamed(rule, "");
+    fRules->Add(name);
+}
+
+void MHMatrix::AddColumns(const MHMatrix *matrix)
+{
+    TIter Next(matrix->fRules);
+    TObject *rule=NULL;
+    while ((rule=Next()))
+        AddColumn(rule->GetName());
 }
 
@@ -264 +283 @@
 {
     int n=0;
-    TIter Next(fData);
+
+    TIter Next(fData->GetSize() ? fData : fRules);
     MData *data = NULL;
     while ((data=(MData*)Next()))
     {
-        *fLog << all << " Column " << setw(3) << n++ << ": ";
+        *fLog << all << " Column " << setw(3) << n++ << ": " << flush;
         data->Print();
         *fLog << endl;
@@ -275 +295 @@
     fM.Print();
 }
+
+const TMatrix *MHMatrix::InvertPosDef()
+{
+    /*
+     ----------------------------------
+      Substract Mean of Rows from Rows
+     ----------------------------------
+
+    const Int_t rows = fM.GetNrows();
+    const Int_t cols = fM.GetNcols();
+
+    for (int i=0; i<rows; i++)
+    {
+        Double_t mean = 0;
+        for (int j=0; j<cols; j++)
+            mean += fM(i, j);
+        mean /= cols;
+
+        for (int j=0; j<cols; j++)
+            fM(i, j) -= mean;
+    }
+    */
+    /*
+     ----------------------------------
+      Substract Mean of Cols from Cols
+     ----------------------------------
+
+    const Int_t rows = fM.GetNrows();
+    const Int_t cols = fM.GetNcols();
+
+    for (int i=0; i<cols; i++)
+    {
+        Double_t mean = 0;
+        for (int j=0; j<rows; j++)
+            mean += fM(j, i);
+        mean /= rows;
+
+        for (int j=0; j<rows; j++)
+            fM(j, i) -= mean;
+    }
+    */
+
+    TMatrix *m2 = new TMatrix(fM, TMatrix::kTransposeMult, fM);
+
+    Double_t det;
+    m2->Invert(&det);
+    if (det==0)
+    {
+        *fLog << err << "ERROR - MHMatrix::InvertPosDef failed (Matrix is sigular)." << endl;
+        delete m2;
+        return NULL;
+    }
+
+    // m2->Print();
+
+    return m2;
+}
+
+Double_t MHMatrix::CalcDist(const TMatrix &m, const TVector &evt, Int_t num) const
+{
+    const Int_t rows = fM.GetNrows();
+    const Int_t cols = fM.GetNcols();
+
+    TArrayD dists(rows);
+
+    //
+    // Calculate:  v^T * M * v
+    //
+    for (int i=0; i<rows; i++)
+    {
+        TVector col(cols);
+        col = TMatrixRow(fM, i);
+
+        TVector d = evt;
+        d -= col;
+
+        TVector d2 = d;
+        d2 *= m;
+
+        dists[i] = fabs(d2*d);
+    }
+
+    TArrayI idx(rows);
+    TMath::Sort(dists.GetSize(), dists.GetArray(), idx.GetArray(), kFALSE);
+
+    const Int_t n = num<rows ? num : rows;
+
+    Double_t res = 0;
+    for (int i=0; i<n; i++)
+        res += dists[idx[i]];
+
+    return res/n;
+}
+
+Double_t MHMatrix::CalcDist(const TVector &evt, Int_t num)
+{
+    if (!fM2.IsValid())
+    {
+        const TMatrix &m = *InvertPosDef();
+        fM2.ResizeTo(m);
+        fM2 = m;
+        delete &m;
+    }
+
+    return CalcDist(fM2, evt, num);
+}
+
+void MHMatrix::Reassign()
+{
+    TMatrix m = fM;
+    fM.ResizeTo(1,1);
+    fM.ResizeTo(m);
+    fM = m;
+
+    TIter Next(fRules);
+    TObject *obj = NULL;
+    while ((obj=Next()))
+        fData->Add(new MDataChain(obj->GetName()));
+}

trunk/MagicSoft/Mars/mhist/MHMatrix.h

@@ -14 +14 @@
 {
 protected:
-    Int_t   fNumRow; // Number of dimensions of histogram
+    Int_t   fNumRow; //! Number of dimensions of histogram
     TMatrix fM;      // Matrix to be filled
 
-    TList  *fData;   // List of data members (columns)
+    TMatrix fM2;     //!
+
+    TList  *fData;   //! List of data members (columns)
+    TList  *fRules;  // List of data members as text for storage
 
     void AddRow();
@@ -26 +29 @@
 
 public:
-    MHMatrix();
+    MHMatrix(const char *name=NULL, const char *title=NULL);
     ~MHMatrix();
 
     void AddColumn(const char *name);
+    void AddColumns(const MHMatrix *matrix);
 
-    TMatrix &GetM() { return fM; }
+    //    TMatrix &GetM() { return fM; }
     const TMatrix &GetM() const { return fM; }
+    const TList *GetRules() const { return fRules; }
 
     //void Draw(Option_t *opt=NULL);
@@ -39 +44 @@
     void Print(Option_t *) const;
 
+    const TMatrix *InvertPosDef();
+
+    Double_t CalcDist(const TMatrix &m, const TVector &v, Int_t num = 25) const;
+    Double_t CalcDist(const TVector &v, Int_t num = 25);
+
+    void Reassign();
+
     ClassDef(MHMatrix, 1) // Multidimensional Matrix to store events
 };

trunk/MagicSoft/Mars/mhist/Makefile

@@ -48 +48 @@
            MHTimeDiffTime.cc \
            MHTimeDiffTheta.cc \
+           MHCompProb.cc \
+           MHHadroness.cc \
            MHMcEnergy.cc \
            MHMcEfficiency.cc \
@@ -53 +55 @@
            MHMcEfficiencyEnergy.cc \
            MHMcEnergyImpact.cc \
+           MHMcEnergyMigration.cc \
            MHThetabarTime.cc \
            MHMcRate.cc \