source: trunk/Mars/manalysis/MMultiDimDistCalc.cc@ 13449

Last change on this file since 13449 was 8907, checked in by tbretz, 16 years ago
*** empty log message ***
File size: 6.7 KB
Line 
1/* ======================================================================== *\
2!
3! *
4! * This file is part of MARS, the MAGIC Analysis and Reconstruction
5! * Software. It is distributed to you in the hope that it can be a useful
6! * and timesaving tool in analysing Data of imaging Cerenkov telescopes.
7! * It is distributed WITHOUT ANY WARRANTY.
8! *
9! * Permission to use, copy, modify and distribute this software and its
10! * documentation for any purpose is hereby granted without fee,
11! * provided that the above copyright notice appear in all copies and
12! * that both that copyright notice and this permission notice appear
13! * in supporting documentation. It is provided "as is" without express
14! * or implied warranty.
15! *
16!
17!
18! Author(s): Thomas Bretz, 5/2002 <mailto:tbretz@astro.uni-wuerzburg.de>
19! Author(s): Rudy Bock, 5/2002 <mailto:rkb@mppmu.mpg.de>
20!
21! Copyright: MAGIC Software Development, 2000-2002
22!
23!
24\* ======================================================================== */
25
26/////////////////////////////////////////////////////////////////////////////
27//
28// MMultiDimDistCalc
29//
30// Calculated a multidimensional distance. It calculates the distance to
31// all vectors in a given matrix describing Gammas and another one
32// describing Hadrons (non gammas). The shortest distances are avaraged.
33// How many distances are used for avaraging can be specified in the
34// constructor.
35//
36// * If you want to use the nearest neighbor function for calculation use:
37// MMultiDimDistCalc::SetUseKernelMethod(kFALSE);
38// * If you want to use the kernel function for calculation use:
39// MMultiDimDistCalc::SetUseKernelMethod(kTRUE); <default>
40// * To use only the n next neighbors for your calculation use:
41// MMultiDimDistCalc::SetUseNumRows(n);
42// * To use all reference events set the number to 0 <default>
43//
44////////////////////////////////////////////////////////////////////////////
45#include "MMultiDimDistCalc.h"
46
47#include <fstream>
48
49#include <math.h>
50
51#include <TVector.h>
52
53#include "MHMatrix.h" // must be before MLogManip.h
54
55#include "MLog.h"
56#include "MLogManip.h"
57
58#include "MParList.h"
59//#include "MDataPhrase.h"
60#include "MDataArray.h"
61
62#include "MParameters.h"
63
64ClassImp(MMultiDimDistCalc);
65
66using namespace std;
67
68static const TString gsDefName = "MMultiDimDistCalc";
69static const TString gsDefTitle = "Calculate Hadronness with Nearest Neighbor/Kernel Method";
70
71// --------------------------------------------------------------------------
72//
73// Setup histograms and the number of distances which are used for
74// avaraging in CalcDist
75//
76MMultiDimDistCalc::MMultiDimDistCalc(const char *name, const char *title)
77 : fNum(0), fUseKernel(kTRUE), fHadronnessName("MHadronness"), fData(NULL)
78{
79 //
80 // set the name and title of this object
81 //
82 fName = name ? name : gsDefName.Data();
83 fTitle = title ? title : gsDefTitle.Data();
84
85 /*
86 fData = new TList;
87 fData->SetOwner();
88 */
89}
90
91// --------------------------------------------------------------------------
92//
93// Needs:
94// - MatrixGammas [MHMatrix]
95// - MatrixHadrons {MHMatrix]
96// - MHadroness
97// - all data containers used to build the matrixes
98//
99// The matrix object can be filles using MFillH. And must be of the same
100// number of columns (with the same meaning).
101//
102Int_t MMultiDimDistCalc::PreProcess(MParList *plist)
103{
104 fMGammas = (MHMatrix*)plist->FindObject("MatrixGammas", "MHMatrix");
105 if (!fMGammas)
106 {
107 *fLog << err << dbginf << "MatrixGammas [MHMatrix] not found... aborting." << endl;
108 return kFALSE;
109 }
110
111 fMHadrons = (MHMatrix*)plist->FindObject("MatrixHadrons", "MHMatrix");
112 if (!fMHadrons)
113 {
114 *fLog << err << dbginf << "MatrixHadrons [MHMatrix] not found... aborting." << endl;
115 return kFALSE;
116 }
117
118 if (fMGammas->GetM().GetNcols() != fMHadrons->GetM().GetNcols())
119 {
120 *fLog << err << dbginf << "Error matrices have different numbers of columns... aborting." << endl;
121 return kFALSE;
122 }
123
124 /*
125 TIter Next(fMGammas->GetRules());
126 TObject *data=NULL;
127 while ((data=Next()))
128 {
129 MDataPhrase *chain = new MDataPhrase(data->GetName());
130 if (!chain->PreProcess(plist))
131 {
132 delete chain;
133 return kFALSE;
134 }
135 fData->Add(chain);
136 }
137 */
138 fData = fMGammas->GetColumns();
139 if (!fData)
140 {
141 *fLog << err << dbginf << "Error matrix doesn't contain columns... aborting." << endl;
142 return kFALSE;
143 }
144
145 if (!fData->PreProcess(plist))
146 {
147 *fLog << err << dbginf << "PreProcessing of the MDataArray failed for the columns failed... aborting." << endl;
148 return kFALSE;
149 }
150
151 fHadronness = (MParameterD*)plist->FindCreateObj("MParameterD", fHadronnessName);
152 if (!fHadronness)
153 return kFALSE;
154
155 return kTRUE;
156}
157
158// --------------------------------------------------------------------------
159//
160// Evaluates the avarage over the fNum shortest distances in a
161// multidimensional space to a matrix (set of vectors) describing
162// gammas and hadrons.
163// The hadroness of the event is defines as the avarage distance to the
164// set of gammas (dg) divided by the avarage distance to the set of
165// hadrons (dh). Because this value is in teh range [0, inf] it is
166// transformed into [0,1] by:
167// H = exp(-dh/dg);
168//
169Int_t MMultiDimDistCalc::Process()
170{
171 // first copy the data from the data array to a vector event
172 TVector event;
173 *fData >> event;
174
175 Int_t numg = fNum;
176 Int_t numh = fNum;
177 if (fNum==0)
178 {
179 numg = fMGammas->GetM().GetNrows();
180 numh = fMHadrons->GetM().GetNrows();
181 }
182 if (fUseKernel)
183 {
184 numg = -numg;
185 numh = -numh;
186 }
187
188 const Double_t dg = fMGammas->CalcDist(event, numg);
189 const Double_t dh = fMHadrons->CalcDist(event, numh);
190
191 if (dg<0 || dh<0)
192 {
193 *fLog << err << "MHMatrix::CalcDist failed (dg=" << dg << ", dh=" << dh << ")... aborting" << endl;
194 return kERROR;
195 }
196
197 fHadronness->SetVal(dg==0 ? 0 : exp(-dh/dg));
198 fHadronness->SetReadyToSave();
199
200 return kTRUE;
201}
202
203void MMultiDimDistCalc::StreamPrimitive(ostream &out) const
204{
205 out << " MMultiDimDist " << GetUniqueName();
206
207 if (fName!=gsDefName || fTitle!=gsDefTitle)
208 {
209 out << "(\"" << fName << "\"";
210 if (fTitle!=gsDefTitle)
211 out << ", \"" << fTitle << "\")";
212 }
213 out << ";" << endl;
214
215 if (fHadronnessName!="MHadronness")
216 out << " " << GetUniqueName() << ".SetHadronnessName(\"" << fHadronnessName << "\");" << endl;
217 if (fNum!=0)
218 out << " " << GetUniqueName() << ".SetUseNumRows(" << fNum << ");" << endl;
219 if (fUseKernel)
220 out << " " << GetUniqueName() << ".SetUseKernelMethod();" << endl;
221}
222
Note: See TracBrowser for help on using the repository browser.