]>
Commit | Line | Data |
---|---|---|
81f1d3ae | 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 | ||
4d610fd5 | 18 | //----------------------------------------------------------------------------- |
19 | /// \class AliMUONAlignmentTask | |
20 | /// AliAnalysisTask to align the MUON spectrometer. | |
21 | /// The Task reads as input ESDs and feeds the MUONTracks to AliMUONAlignment. | |
22 | /// The alignment itself is performed by AliMillepede. | |
23 | /// A OCDB entry is written with the alignment parameters. | |
24 | /// A list of graph are written to monitor the alignment parameters. | |
25 | /// | |
26 | /// \author Javier Castillo, CEA/Saclay - Irfu/SPhN | |
27 | //----------------------------------------------------------------------------- | |
28 | ||
29 | #include <fstream> | |
30 | ||
31 | #include <TString.h> | |
32 | #include <TError.h> | |
33 | #include <TGraphErrors.h> | |
34 | #include <TTree.h> | |
35 | #include <TChain.h> | |
36 | #include <TClonesArray.h> | |
37 | #include <TGeoGlobalMagField.h> | |
38 | #include <Riostream.h> | |
39 | ||
40 | #include "AliAnalysisTask.h" | |
41 | #include "AliAnalysisManager.h" | |
42 | #include "AliESDInputHandler.h" | |
43 | #include "AliESDEvent.h" | |
44 | #include "AliESDMuonTrack.h" | |
45 | #include "AliMagF.h" | |
46 | #include "AliCDBManager.h" | |
47 | #include "AliCDBMetaData.h" | |
48 | #include "AliCDBId.h" | |
49 | #include "AliGeomManager.h" | |
50 | ||
51 | #include "AliMpCDB.h" | |
52 | #include "AliMUONAlignment.h" | |
53 | #include "AliMUONTrack.h" | |
54 | #include "AliMUONTrackExtrap.h" | |
55 | #include "AliMUONTrackParam.h" | |
56 | #include "AliMUONGeometryTransformer.h" | |
57 | #include "AliMUONESDInterface.h" | |
58 | ||
59 | #include "AliMUONAlignmentTask.h" | |
60 | ||
61 | ///\cond CLASSIMP | |
62 | ClassImp(AliMUONAlignmentTask) | |
63 | ///\endcond | |
64 | ||
65 | // //________________________________________________________________________ | |
66 | // AliMUONAlignmentTask::AliMUONAlignmentTask(const char *name) | |
67 | // : AliAnalysisTask(name, ""), | |
68 | // fESD(0x0), | |
69 | // fAlign(0x0), | |
70 | // fGeoFilename(""), | |
71 | // fTransform(0x0), | |
72 | // fTrackTot(0), | |
73 | // fTrackOk(0), | |
74 | // fMSDEx(0x0), | |
75 | // fMSDEy(0x0), | |
76 | // fMSDEz(0x0), | |
77 | // fMSDEp(0x0), | |
78 | // fList(0x0) | |
79 | // { | |
80 | // /// Default Constructor | |
81 | // // Define input and output slots here | |
82 | // // Input slot #0 works with a TChain | |
83 | // DefineInput(0, TChain::Class()); | |
84 | // // Output slot #0 writes NTuple/histos into a TList | |
85 | // DefineOutput(0, TList::Class()); | |
86 | ||
87 | // // initialize parameters ... | |
88 | // for(Int_t k=0;k<4*156;k++) { | |
89 | // fParameters[k]=0.; | |
90 | // fErrors[k]=0.; | |
91 | // fPulls[k]=0.; | |
92 | // } | |
93 | ||
94 | // fAlign = new AliMUONAlignment(); | |
95 | // fTransform = new AliMUONGeometryTransformer(); | |
96 | // } | |
97 | ||
98 | //________________________________________________________________________ | |
99 | AliMUONAlignmentTask::AliMUONAlignmentTask(const char *name, const char *geofilename) | |
100 | : AliAnalysisTask(name, ""), | |
101 | fESD(0x0), | |
102 | fAlign(0x0), | |
103 | fGeoFilename(geofilename), | |
104 | fTransform(0x0), | |
105 | fTrackTot(0), | |
106 | fTrackOk(0), | |
107 | fMSDEx(0x0), | |
108 | fMSDEy(0x0), | |
109 | fMSDEz(0x0), | |
110 | fMSDEp(0x0), | |
111 | fList(0x0) | |
112 | { | |
113 | /// Default Constructor | |
114 | // Define input and output slots here | |
115 | // Input slot #0 works with a TChain | |
116 | DefineInput(0, TChain::Class()); | |
117 | // Output slot #0 writes NTuple/histos into a TList | |
118 | DefineOutput(0, TList::Class()); | |
119 | ||
120 | // initialize parameters ... | |
121 | for(Int_t k=0;k<4*156;k++) { | |
122 | fParameters[k]=0.; | |
123 | fErrors[k]=0.; | |
124 | fPulls[k]=0.; | |
125 | } | |
126 | ||
127 | fAlign = new AliMUONAlignment(); | |
128 | fTransform = new AliMUONGeometryTransformer(); | |
129 | } | |
130 | ||
131 | ||
132 | //________________________________________________________________________ | |
133 | AliMUONAlignmentTask::AliMUONAlignmentTask(const AliMUONAlignmentTask& obj) | |
134 | : AliAnalysisTask(obj), | |
135 | fESD(0x0), | |
136 | fAlign(0x0), | |
137 | fGeoFilename(""), | |
138 | fTransform(0x0), | |
139 | fTrackTot(0), | |
140 | fTrackOk(0), | |
141 | fMSDEx(0x0), | |
142 | fMSDEy(0x0), | |
143 | fMSDEz(0x0), | |
144 | fMSDEp(0x0), | |
145 | fList(0x0) | |
146 | { | |
147 | /// Copy constructor | |
148 | fESD = obj.fESD; | |
149 | fAlign = obj.fAlign; | |
150 | fGeoFilename = obj.fGeoFilename; | |
151 | fTransform = obj.fTransform; | |
152 | fTrackTot = obj.fTrackTot; | |
153 | fTrackOk = obj.fTrackOk; | |
154 | fMSDEx = obj.fMSDEx; | |
155 | fMSDEy = obj.fMSDEy; | |
156 | fMSDEz = obj.fMSDEz; | |
157 | fMSDEp = obj.fMSDEp; | |
158 | fList = obj.fList; | |
159 | } | |
160 | ||
161 | //________________________________________________________________________ | |
162 | AliMUONAlignmentTask& AliMUONAlignmentTask::operator=(const AliMUONAlignmentTask& other) | |
163 | { | |
164 | /// Assignment | |
165 | AliAnalysisTask::operator=(other); | |
166 | fESD = other.fESD; | |
167 | fAlign = other.fAlign; | |
168 | fGeoFilename = other.fGeoFilename; | |
169 | fTransform = other.fTransform; | |
170 | fTrackTot = other.fTrackTot; | |
171 | fTrackOk = other.fTrackOk; | |
172 | fMSDEx = other.fMSDEx; | |
173 | fMSDEy = other.fMSDEy; | |
174 | fMSDEz = other.fMSDEz; | |
175 | fMSDEp = other.fMSDEp; | |
176 | fList = other.fList; | |
177 | ||
178 | return *this; | |
179 | } | |
180 | ||
181 | //________________________________________________________________________ | |
182 | AliMUONAlignmentTask::~AliMUONAlignmentTask() | |
183 | { | |
184 | /// Destructor | |
185 | if (fAlign) delete fAlign; | |
186 | if (fTransform) delete fTransform; | |
187 | } | |
188 | ||
189 | //________________________________________________________________________ | |
190 | void AliMUONAlignmentTask::LocalInit() | |
191 | { | |
192 | /// Local initialization, called once per task on the client machine | |
193 | /// where the analysis train is assembled | |
194 | AliMpCDB::LoadMpSegmentation(); | |
195 | ||
196 | // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex) | |
197 | if ( ! AliGeomManager::GetGeometry() ) { | |
198 | AliGeomManager::LoadGeometry(fGeoFilename.Data()); | |
199 | if (! AliGeomManager::GetGeometry() ) { | |
200 | Error("MUONAlignment", "getting geometry from file %s failed", fGeoFilename.Data()); | |
201 | return; | |
202 | } | |
203 | } | |
204 | ||
205 | // set mag field | |
206 | // waiting for mag field in CDB | |
207 | if (!TGeoGlobalMagField::Instance()->GetField()) { | |
208 | printf("Loading field map...\n"); | |
209 | AliMagF* field = new AliMagF("Maps","Maps",2,0.,0., 10.,AliMagF::k5kG); | |
210 | TGeoGlobalMagField::Instance()->SetField(field); | |
211 | } | |
212 | // set the magnetic field for track extrapolations | |
213 | AliMUONTrackExtrap::SetField(); | |
214 | ||
215 | // Set initial values here, good guess may help convergence | |
216 | // St 1 | |
217 | // Int_t iPar = 0; | |
218 | // fParameters[iPar++] = 0.010300 ; fParameters[iPar++] = 0.010600 ; fParameters[iPar++] = 0.000396 ; | |
219 | ||
220 | ||
221 | fAlign->InitGlobalParameters(fParameters); | |
222 | ||
223 | ||
224 | fTransform->LoadGeometryData(); | |
225 | fAlign->SetGeometryTransformer(fTransform); | |
226 | ||
227 | // Do alignment with magnetic field off | |
228 | fAlign->SetBFieldOn(kFALSE); | |
229 | ||
230 | // Set tracking station to use | |
231 | // Bool_t bStOnOff[5] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE}; | |
232 | Bool_t bChOnOff[10] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kFALSE,kTRUE,kTRUE,kTRUE,kTRUE}; | |
233 | ||
234 | // Set degrees of freedom to align (see AliMUONAlignment) | |
235 | fAlign->AllowVariations(bChOnOff); | |
236 | ||
237 | // Fix parameters or add constraints here | |
238 | // for (Int_t iSt=0; iSt<5; iSt++) | |
239 | // if (!bStOnOff[iSt]) fAlign->FixStation(iSt+1); | |
240 | for (Int_t iCh=0; iCh<10; iCh++) | |
241 | if (!bChOnOff[iCh]) fAlign->FixChamber(iCh+1); | |
242 | ||
243 | // Left and right sides of the detector are independent, one can choose to align | |
244 | // only one side | |
245 | Bool_t bSpecLROnOff[2] = {kTRUE,kTRUE}; | |
246 | fAlign->FixHalfSpectrometer(bChOnOff,bSpecLROnOff); | |
247 | ||
248 | fAlign->SetChOnOff(bChOnOff); | |
249 | fAlign->SetSpecLROnOff(bChOnOff); | |
250 | ||
251 | // Here we can fix some detection elements | |
252 | fAlign->FixDetElem(908); | |
253 | fAlign->FixDetElem(1020); | |
254 | ||
255 | // Set predifined global constrains: X, Y, P, XvsZ, YvsZ, PvsZ, XvsY, YvsY, PvsY | |
256 | // Bool_t bVarXYT[9] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE}; | |
257 | // Bool_t bDetTLBR[4] = {kFALSE,kTRUE,kFALSE,kTRUE}; | |
258 | // fAlign->AddConstraints(bChOnOff,bVarXYT,bDetTLBR,bSpecLROnOff); | |
259 | ||
260 | } | |
261 | ||
262 | //________________________________________________________________________ | |
263 | void AliMUONAlignmentTask::ConnectInputData(Option_t *) | |
264 | { | |
265 | /// Connect ESD here. Called on each input data change. | |
266 | ||
267 | // Connect ESD here | |
268 | TTree* esdTree = dynamic_cast<TTree*> (GetInputData(0)); | |
269 | if (!esdTree) { | |
270 | Printf("ERROR: Could not read chain from input slot 0"); | |
271 | } | |
272 | else { | |
273 | AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()); | |
274 | if (!esdH) { | |
275 | Printf("ERROR: Could not get ESDInputHandler"); | |
276 | } | |
277 | else { | |
278 | fESD = esdH->GetEvent(); | |
279 | } | |
280 | } | |
281 | } | |
282 | ||
283 | //________________________________________________________________________ | |
284 | void AliMUONAlignmentTask::CreateOutputObjects() | |
285 | { | |
286 | /// Executed once on each worker (machine actually running the analysis code) | |
287 | // | |
288 | // This method has to be called INSIDE the user redefined CreateOutputObjects | |
289 | // method, before creating each object corresponding to the output containers | |
290 | // that are to be written to a file. This need to be done in general for the big output | |
291 | // objects that may not fit memory during processing. | |
292 | // OpenFile(0); | |
293 | ||
294 | // Creating graphs | |
295 | fMSDEx = new TGraphErrors(156); | |
296 | fMSDEy = new TGraphErrors(156); | |
297 | fMSDEz = new TGraphErrors(156); | |
298 | fMSDEp = new TGraphErrors(156); | |
299 | ||
300 | // Add Ntuples to the list | |
301 | fList = new TList(); | |
302 | fList->Add(fMSDEx); | |
303 | fList->Add(fMSDEy); | |
304 | fList->Add(fMSDEz); | |
305 | fList->Add(fMSDEp); | |
306 | } | |
307 | ||
308 | //________________________________________________________________________ | |
309 | void AliMUONAlignmentTask::Exec(Option_t *) | |
310 | { | |
311 | /// Main loop, called for each event | |
312 | if (!fESD) { | |
313 | Printf("ERROR: fESD not available"); | |
314 | return; | |
315 | } | |
316 | ||
317 | Double_t trackParams[8] = {0.,0.,0.,0.,0.,0.,0.,0.}; | |
318 | ||
319 | ||
320 | ||
321 | Int_t nTracks = Int_t(fESD->GetNumberOfMuonTracks()); | |
322 | // if (!event%100) cout << " there are " << nTracks << " tracks in event " << event << endl; | |
323 | for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) { | |
324 | AliESDMuonTrack* esdTrack = fESD->GetMuonTrack(iTrack); | |
325 | if (!esdTrack->ClustersStored()) continue; | |
326 | Double_t invBenMom = esdTrack->GetInverseBendingMomentum(); | |
327 | // fInvBenMom->Fill(invBenMom); | |
328 | // fBenMom->Fill(1./invBenMom); | |
329 | if (TMath::Abs(invBenMom)<=1.04) { | |
330 | AliMUONTrack track; | |
331 | AliMUONESDInterface::ESDToMUON(*esdTrack, track); | |
332 | fAlign->ProcessTrack(&track); | |
333 | fAlign->LocalFit(fTrackOk++,trackParams,0); | |
334 | } | |
335 | fTrackTot++; | |
336 | } | |
337 | ||
338 | // Post final data. Write histo list to a file with option "RECREATE" | |
339 | PostData(0,fList); | |
340 | } | |
341 | ||
342 | //________________________________________________________________________ | |
343 | void AliMUONAlignmentTask::Terminate(const Option_t*) | |
344 | { | |
345 | /// Called once per task on the client machine at the end of the analysis. | |
346 | ||
347 | cout << "Processed " << fTrackTot << " Tracks." << endl; | |
348 | // Perform global fit | |
349 | fAlign->GlobalFit(fParameters,fErrors,fPulls); | |
350 | ||
351 | cout << "Done with GlobalFit " << endl; | |
352 | ||
353 | // Update pointers reading them from the output slot | |
354 | fList = (TList*)GetOutputData(0); | |
355 | fMSDEx = (TGraphErrors*)fList->At(0); | |
356 | fMSDEy = (TGraphErrors*)fList->At(1); | |
357 | fMSDEz = (TGraphErrors*)fList->At(2); | |
358 | fMSDEp = (TGraphErrors*)fList->At(3); | |
359 | ||
360 | // Store results | |
361 | Double_t DEid[156] = {0}; | |
362 | Double_t MSDEx[156] = {0}; | |
363 | Double_t MSDEy[156] = {0}; | |
364 | Double_t MSDEz[156] = {0}; | |
365 | Double_t MSDEp[156] = {0}; | |
366 | Double_t DEidErr[156] = {0}; | |
367 | Double_t MSDExErr[156] = {0}; | |
368 | Double_t MSDEyErr[156] = {0}; | |
369 | Double_t MSDEzErr[156] = {0}; | |
370 | Double_t MSDEpErr[156] = {0}; | |
371 | Int_t lNDetElem = 4*2+4*2+18*2+26*2+26*2; | |
372 | Int_t lNDetElemCh[10] = {4,4,4,4,18,18,26,26,26,26}; | |
373 | // Int_t lSNDetElemCh[10] = {4,8,12,16,34,52,78,104,130,156}; | |
374 | Int_t idOffset = 0; // 400 | |
375 | Int_t lSDetElemCh = 0; | |
376 | for(Int_t iDE=0; iDE<lNDetElem; iDE++){ | |
377 | DEidErr[iDE] = 0.; | |
378 | DEid[iDE] = idOffset+100; | |
379 | DEid[iDE] += iDE; | |
380 | lSDetElemCh = 0; | |
381 | for(Int_t iCh=0; iCh<9; iCh++){ | |
382 | lSDetElemCh += lNDetElemCh[iCh]; | |
383 | if (iDE>=lSDetElemCh) { | |
384 | DEid[iDE] += 100; | |
385 | DEid[iDE] -= lNDetElemCh[iCh]; | |
386 | } | |
387 | } | |
388 | MSDEx[iDE]=fParameters[3*iDE+0]; | |
389 | MSDEy[iDE]=fParameters[3*iDE+1]; | |
390 | MSDEz[iDE]=fParameters[3*iDE+3]; | |
391 | MSDEp[iDE]=fParameters[3*iDE+2]; | |
392 | MSDExErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+0); | |
393 | MSDEyErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+1); | |
394 | MSDEzErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+3); | |
395 | MSDEpErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+2); | |
396 | fMSDEx->SetPoint(iDE,DEid[iDE],fParameters[3*iDE+0]); | |
397 | fMSDEx->SetPoint(iDE,DEidErr[iDE],(Double_t)fAlign->GetParError(3*iDE+0)); | |
398 | fMSDEy->SetPoint(iDE,DEid[iDE],fParameters[3*iDE+1]); | |
399 | fMSDEy->SetPoint(iDE,DEidErr[iDE],(Double_t)fAlign->GetParError(3*iDE+1)); | |
400 | fMSDEp->SetPoint(iDE,DEid[iDE],fParameters[3*iDE+2]); | |
401 | fMSDEp->SetPoint(iDE,DEidErr[iDE],(Double_t)fAlign->GetParError(3*iDE+2)); | |
402 | fMSDEz->SetPoint(iDE,DEid[iDE],fParameters[3*iDE+3]); | |
403 | fMSDEz->SetPoint(iDE,DEidErr[iDE],(Double_t)fAlign->GetParError(3*iDE+3)); | |
404 | } | |
405 | ||
406 | // Post final data. Write histo list to a file with option "RECREATE" | |
407 | PostData(0,fList); | |
408 | ||
409 | ||
410 | // Re Align | |
411 | AliMUONGeometryTransformer *newTransform = fAlign->ReAlign(fTransform,fParameters,true); | |
412 | newTransform->WriteTransformations("transform2ReAlign.dat"); | |
413 | ||
414 | // Generate realigned data in local cdb | |
415 | const TClonesArray* array = newTransform->GetMisAlignmentData(); | |
416 | ||
417 | // 100 mum residual resolution for chamber misalignments? | |
418 | fAlign->SetAlignmentResolution(array,-1,0.01,0.01,0.004,0.003); | |
419 | ||
420 | // CDB manager | |
421 | AliCDBManager* cdbManager = AliCDBManager::Instance(); | |
422 | cdbManager->SetDefaultStorage("local://ReAlignCDB"); | |
423 | ||
424 | AliCDBMetaData* cdbData = new AliCDBMetaData(); | |
425 | cdbData->SetResponsible("Dimuon Offline project"); | |
426 | cdbData->SetComment("MUON alignment objects with residual misalignment"); | |
427 | AliCDBId id("MUON/Align/Data", 0, AliCDBRunRange::Infinity()); | |
428 | cdbManager->Put(const_cast<TClonesArray*>(array), id, cdbData); | |
429 | ||
430 | } | |
431 |