Corrected Doxygen warnings:
[u/mrichter/AliRoot.git] / MUON / AliMUONClusterFinderSimpleFit.cxx
CommitLineData
e036a7af 1/**************************************************************************
2* Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3* *
4* Author: The ALICE Off-line Project. *
5* Contributors are mentioned in the code where appropriate. *
6* *
7* Permission to use, copy, modify and distribute this software and its *
8* documentation strictly for non-commercial purposes is hereby granted *
9* without fee, provided that the above copyright notice appears in all *
10* copies and that both the copyright notice and this permission notice *
11* appear in the supporting documentation. The authors make no claims *
12* about the suitability of this software for any purpose. It is *
13* provided "as is" without express or implied warranty. *
14**************************************************************************/
15
16// $Id$
17
18#include "AliMUONClusterFinderSimpleFit.h"
19
20#include "AliLog.h"
21#include "AliMpDEManager.h"
22#include "AliMpStationType.h"
23#include "AliMUONCluster.h"
24#include "AliMUONConstants.h"
25#include "AliMUONDigit.h"
26#include "AliMUONMathieson.h"
27#include "AliMUONPad.h"
28#include "AliMUONClusterFinderCOG.h"
29#include "AliMpArea.h"
30#include "TClonesArray.h"
31#include "TObjArray.h"
32#include "TVector2.h"
33#include "TVirtualFitter.h"
34#include "TF1.h"
35
36/// \class AliMUONClusterFinderSimpleFit
37///
38/// Basic cluster finder
39///
40/// We simply use AliMUONPreClusterFinder to get basic cluster,
41/// and then we try to fit the charge repartition using a Mathieson
42/// distribution, varying the position.
43///
44/// FIXME: this one is still at the developping stage...
45///
46/// \author Laurent Aphecetche
47
78649106 48/// \cond CLASSIMP
e036a7af 49ClassImp(AliMUONClusterFinderSimpleFit)
50/// \endcond
51
52namespace
53{
54 //___________________________________________________________________________
55 void
56 FitFunction(Int_t& /*notused*/, Double_t* /*notused*/,
57 Double_t& f, Double_t* par,
58 Int_t /*notused*/)
59 {
60 /// Chi2 Function to minimize: Mathieson charge distribution in 2 dimensions
61
62 TObjArray* userObjects = static_cast<TObjArray*>(TVirtualFitter::GetFitter()->GetObjectFit());
63
64 AliMUONCluster* cluster = static_cast<AliMUONCluster*>(userObjects->At(0));
65 AliMUONMathieson* mathieson = static_cast<AliMUONMathieson*>(userObjects->At(1));
66
67 f = 0.0;
68 Float_t qTot = cluster->Charge();
32074b02 69// Float_t chargeCorrel[] = { cluster->Charge(0)/qTot, cluster->Charge(1)/qTot };
70// Float_t qRatio[] = { 1.0/par[2], par[2] };
e036a7af 71
72 for ( Int_t i = 0 ; i < cluster->Multiplicity(); ++i )
73 {
74 AliMUONPad* pad = cluster->Pad(i);
75 // skip pads w/ saturation or other problem(s)
76 if ( pad->Status() ) continue;
ff331318 77 TVector2 lowerLeft = TVector2(par[0],par[1]) - pad->Position() - pad->Dimensions();
e036a7af 78 TVector2 upperRight(lowerLeft + pad->Dimensions()*2.0);
32074b02 79 Float_t estimatedCharge = mathieson->IntXY(lowerLeft.X(),lowerLeft.Y(),
80 upperRight.X(),upperRight.Y());
81// estimatedCharge *= 2/(1+qRatio[pad->Cathode()]);
82 Float_t actualCharge = pad->Charge()/qTot;
e036a7af 83
84 Float_t delta = (estimatedCharge - actualCharge);
85
32074b02 86 f += delta*delta;
e036a7af 87 }
88 }
89}
90
91//_____________________________________________________________________________
92AliMUONClusterFinderSimpleFit::AliMUONClusterFinderSimpleFit()
93: AliMUONVClusterFinder(),
94fClusterFinder(0x0),
95fMathieson(0x0)
96{
97 /// ctor
98}
99
100//_____________________________________________________________________________
101AliMUONClusterFinderSimpleFit::~AliMUONClusterFinderSimpleFit()
102{
103 /// dtor
104 delete fClusterFinder;
105 delete fMathieson;
106}
107
108//_____________________________________________________________________________
109Bool_t
110AliMUONClusterFinderSimpleFit::Prepare(const AliMpVSegmentation* segmentations[2],
111 TClonesArray* digits[2])
112{
113 /// Prepare for clustering
114
115 // FIXME: should we get the Mathieson from elsewhere ?
116
117 // Find out the DetElemId
118 Int_t detElemId(-1);
119
120 for ( Int_t i = 0; i < 2; ++i )
121 {
122 AliMUONDigit* d = static_cast<AliMUONDigit*>(digits[i]->First());
123 if (d)
124 {
125 detElemId = d->DetElemId();
126 break;
127 }
128 }
129
130 if ( detElemId < 0 )
131 {
132 AliWarning("Could not find DE. Probably no digits at all ?");
133 return kFALSE;
134 }
135
866c3232 136 AliMp::StationType stationType = AliMpDEManager::GetStationType(detElemId);
e036a7af 137
138 Float_t kx3 = AliMUONConstants::SqrtKx3();
139 Float_t ky3 = AliMUONConstants::SqrtKy3();
140 Float_t pitch = AliMUONConstants::Pitch();
141
866c3232 142 if ( stationType == AliMp::kStation1 )
e036a7af 143 {
144 kx3 = AliMUONConstants::SqrtKx3St1();
145 ky3 = AliMUONConstants::SqrtKy3St1();
146 pitch = AliMUONConstants::PitchSt1();
147 }
148
149 delete fMathieson;
150 fMathieson = new AliMUONMathieson;
151
152 fMathieson->SetPitch(pitch);
153 fMathieson->SetSqrtKx3AndDeriveKx2Kx4(kx3);
154 fMathieson->SetSqrtKy3AndDeriveKy2Ky4(ky3);
155
156 delete fClusterFinder;
157 fClusterFinder = new AliMUONClusterFinderCOG;
158 return fClusterFinder->Prepare(segmentations,digits);
159}
160
161//_____________________________________________________________________________
162AliMUONCluster*
163AliMUONClusterFinderSimpleFit::NextCluster()
164{
165 /// Returns next cluster
166
167 if ( !fClusterFinder ) return 0x0;
168 AliMUONCluster* cluster = fClusterFinder->NextCluster();
169 if ( cluster )
170 {
171 ComputePosition(*cluster);
172
173 if ( cluster->Charge() < 7 )
174 {
175 // skip that one
176 return NextCluster();
177 }
178 }
179 return cluster;
180}
181
182//_____________________________________________________________________________
183void
184AliMUONClusterFinderSimpleFit::ComputePosition(AliMUONCluster& cluster)
185{
186 /// Compute the position of the given cluster, by fitting a Mathieson
187 /// charge distribution to it
188
189 TVirtualFitter* fitter = TVirtualFitter::Fitter(0,2);
190 fitter->SetFCN(FitFunction);
191
192 if ( cluster.Multiplicity() < 3 ) return;
193
194 // We try a Mathieson fit, starting
195 // with the center-of-gravity estimate as a first guess
196 // for the cluster center.
197
198 Double_t xCOG = cluster.Position().X();
199 Double_t yCOG = cluster.Position().Y();
200
201 Float_t stepX = 0.01; // cm
202 Float_t stepY = 0.01; // cm
203
32074b02 204 Double_t arg(-1); // disable printout
e036a7af 205
32074b02 206 fitter->ExecuteCommand("SET PRINT",&arg,1);
e036a7af 207
208 fitter->SetParameter(0,"cluster X position",xCOG,stepX,0,0);
209 fitter->SetParameter(1,"cluster Y position",yCOG,stepY,0,0);
e036a7af 210
211 TObjArray userObjects;
212
213 userObjects.Add(&cluster);
214 userObjects.Add(fMathieson);
215
e036a7af 216 fitter->SetObjectFit(&userObjects);
217
32074b02 218 Int_t val = fitter->ExecuteCommand("MIGRAD",0,0);
219 AliDebug(1,Form("ExecuteCommand returned value=%d",val));
220 if ( val )
221 {
222 // fit failed. Using COG results, with big errors
223 AliWarning("Fit failed. Using COG results for cluster=");
224 StdoutToAliWarning(cluster.Print());
225 cluster.SetPosition(TVector2(xCOG,yCOG),TVector2(TMath::Abs(xCOG),TMath::Abs(yCOG)));
226 cluster.SetChi2(1E3);
227 }
e036a7af 228
229 Double_t results[] = { fitter->GetParameter(0),
230 fitter->GetParameter(1) };
231
232 Double_t errors[] = { fitter->GetParError(0),
233 fitter->GetParError(1) };
234
235 cluster.SetPosition(TVector2(results[0],results[1]),
236 TVector2(errors[0],errors[1]));
237
32074b02 238 Double_t amin, edm, errdef;
239 Int_t nvpar, nparx;
240
241 fitter->GetStats(amin, edm, errdef, nvpar, nparx);
242
243 Double_t chi2 = amin;
e036a7af 244
32074b02 245 AliDebug(1,Form("Cluster fitted to (x,y)=(%e,%e) (xerr,yerr)=(%e,%e) \n chi2=%e ndf=%d",
e036a7af 246 results[0],results[1],
32074b02 247 errors[0],errors[1],chi2,fitter->GetNumberFreeParameters()));
248 cluster.SetChi2(chi2);
e036a7af 249}
250
251
252