1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
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 **************************************************************************/
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.
26 /// \author Javier Castillo, CEA/Saclay - Irfu/SPhN
27 //-----------------------------------------------------------------------------
33 #include <TGraphErrors.h>
36 #include <TClonesArray.h>
37 #include <TGeoGlobalMagField.h>
38 #include <TGeoManager.h>
39 #include <Riostream.h>
41 #include "AliAnalysisTask.h"
42 #include "AliAnalysisManager.h"
43 #include "AliESDInputHandler.h"
44 #include "AliESDEvent.h"
45 #include "AliESDMuonTrack.h"
47 #include "AliCDBManager.h"
48 #include "AliGRPManager.h"
49 #include "AliCDBMetaData.h"
51 #include "AliGeomManager.h"
54 #include "AliMUONAlignment.h"
55 #include "AliMUONTrack.h"
56 #include "AliMUONTrackExtrap.h"
57 #include "AliMUONTrackParam.h"
58 #include "AliMUONGeometryTransformer.h"
59 #include "AliMUONESDInterface.h"
61 #include "AliMUONAlignmentTask.h"
64 ClassImp(AliMUONAlignmentTask)
67 // //________________________________________________________________________
68 // AliMUONAlignmentTask::AliMUONAlignmentTask(const char *name)
69 // : AliAnalysisTask(name, ""),
82 // /// Default Constructor
83 // // Define input and output slots here
84 // // Input slot #0 works with a TChain
85 // DefineInput(0, TChain::Class());
86 // // Output slot #0 writes NTuple/histos into a TList
87 // DefineOutput(0, TList::Class());
89 // // initialize parameters ...
90 // for(Int_t k=0;k<4*156;k++) {
96 // fAlign = new AliMUONAlignment();
97 // fTransform = new AliMUONGeometryTransformer();
100 //________________________________________________________________________
101 AliMUONAlignmentTask::AliMUONAlignmentTask(const char *name, const char *geofilename, const char *defaultocdb, const char *misalignocdb)
102 : AliAnalysisTask(name, ""),
105 fGeoFilename(geofilename),
106 fMisAlignOCDB(misalignocdb),
107 fDefaultOCDB(defaultocdb),
118 /// Default Constructor
119 // Define input and output slots here
120 // Input slot #0 works with a TChain
121 DefineInput(0, TChain::Class());
122 // Output slot #0 writes NTuple/histos into a TList
123 DefineOutput(0, TList::Class());
125 // initialize parameters ...
126 for(Int_t k=0;k<4*156;k++) {
132 fAlign = new AliMUONAlignment();
133 fTransform = new AliMUONGeometryTransformer();
137 //________________________________________________________________________
138 AliMUONAlignmentTask::AliMUONAlignmentTask(const AliMUONAlignmentTask& obj)
139 : AliAnalysisTask(obj),
158 fGeoFilename = obj.fGeoFilename;
159 fTransform = obj.fTransform;
160 fTrackTot = obj.fTrackTot;
161 fTrackOk = obj.fTrackOk;
162 fLastRunNumber = obj.fLastRunNumber;
169 // initialize parameters ...
170 for(Int_t k=0;k<4*156;k++) {
171 fParameters[k]=obj.fParameters[k];
172 fErrors[k]=obj.fErrors[k];
173 fPulls[k]=obj.fPulls[k];
178 //________________________________________________________________________
179 AliMUONAlignmentTask& AliMUONAlignmentTask::operator=(const AliMUONAlignmentTask& other)
182 AliAnalysisTask::operator=(other);
184 fAlign = other.fAlign;
185 fGeoFilename = other.fGeoFilename;
186 fMisAlignOCDB = other.fMisAlignOCDB;
187 fDefaultOCDB = other.fDefaultOCDB;
188 fTransform = other.fTransform;
189 fTrackTot = other.fTrackTot;
190 fTrackOk = other.fTrackOk;
191 fLastRunNumber = other.fLastRunNumber;
192 fMSDEx = other.fMSDEx;
193 fMSDEy = other.fMSDEy;
194 fMSDEz = other.fMSDEz;
195 fMSDEp = other.fMSDEp;
198 // initialize parameters ...
199 for(Int_t k=0;k<4*156;k++) {
200 fParameters[k]=other.fParameters[k];
201 fErrors[k]=other.fErrors[k];
202 fPulls[k]=other.fPulls[k];
208 //________________________________________________________________________
209 AliMUONAlignmentTask::~AliMUONAlignmentTask()
212 if (fAlign) delete fAlign;
213 if (fTransform) delete fTransform;
216 //________________________________________________________________________
217 void AliMUONAlignmentTask::LocalInit()
219 /// Local initialization, called once per task on the client machine
220 /// where the analysis train is assembled
222 // Prepare(fGeoFilename.Data(),fDefaultOCDB.Data(),fMisAlignOCDB.Data());
223 Prepare(fGeoFilename.Data(),"local://$ALICE_ROOT/OCDB",fMisAlignOCDB.Data());
226 // Set initial values here, good guess may help convergence
229 // fParameters[iPar++] = 0.010300 ; fParameters[iPar++] = 0.010600 ; fParameters[iPar++] = 0.000396 ;
232 fAlign->InitGlobalParameters(fParameters);
235 fTransform->LoadGeometryData();
236 fAlign->SetGeometryTransformer(fTransform);
238 // Do alignment with magnetic field off
239 fAlign->SetBFieldOn(kFALSE);
241 // Set tracking station to use
242 // Bool_t bStOnOff[5] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
243 Bool_t bChOnOff[10] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kFALSE,kTRUE,kTRUE,kTRUE,kTRUE};
245 // Set degrees of freedom to align (see AliMUONAlignment)
246 fAlign->AllowVariations(bChOnOff);
248 // Fix parameters or add constraints here
249 // for (Int_t iSt=0; iSt<5; iSt++)
250 // if (!bStOnOff[iSt]) fAlign->FixStation(iSt+1);
251 for (Int_t iCh=0; iCh<10; iCh++)
252 if (!bChOnOff[iCh]) fAlign->FixChamber(iCh+1);
254 // Left and right sides of the detector are independent, one can choose to align
256 Bool_t bSpecLROnOff[2] = {kTRUE,kTRUE};
257 fAlign->FixHalfSpectrometer(bChOnOff,bSpecLROnOff);
259 fAlign->SetChOnOff(bChOnOff);
260 fAlign->SetSpecLROnOff(bChOnOff);
262 // Here we can fix some detection elements
263 fAlign->FixDetElem(908);
264 fAlign->FixDetElem(1020);
266 // Set predifined global constrains: X, Y, P, XvsZ, YvsZ, PvsZ, XvsY, YvsY, PvsY
267 // Bool_t bVarXYT[9] = {kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE,kTRUE};
268 // Bool_t bDetTLBR[4] = {kFALSE,kTRUE,kFALSE,kTRUE};
269 // fAlign->AddConstraints(bChOnOff,bVarXYT,bDetTLBR,bSpecLROnOff);
273 //________________________________________________________________________
274 void AliMUONAlignmentTask::ConnectInputData(Option_t *)
276 /// Connect ESD here. Called on each input data change.
279 TTree* esdTree = dynamic_cast<TTree*> (GetInputData(0));
281 Printf("ERROR: Could not read chain from input slot 0");
284 AliESDInputHandler *esdH = dynamic_cast<AliESDInputHandler*> (AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler());
286 Printf("ERROR: Could not get ESDInputHandler");
289 fESD = esdH->GetEvent();
294 //________________________________________________________________________
295 void AliMUONAlignmentTask::CreateOutputObjects()
297 /// Executed once on each worker (machine actually running the analysis code)
299 // This method has to be called INSIDE the user redefined CreateOutputObjects
300 // method, before creating each object corresponding to the output containers
301 // that are to be written to a file. This need to be done in general for the big output
302 // objects that may not fit memory during processing.
306 fMSDEx = new TGraphErrors(156);
307 fMSDEy = new TGraphErrors(156);
308 fMSDEz = new TGraphErrors(156);
309 fMSDEp = new TGraphErrors(156);
311 // Add Ntuples to the list
319 //________________________________________________________________________
320 void AliMUONAlignmentTask::Exec(Option_t *)
322 /// Main loop, called for each event
324 Printf("ERROR: fESD not available");
328 Double_t trackParams[8] = {0.,0.,0.,0.,0.,0.,0.,0.};
329 if (fESD->GetRunNumber()!=fLastRunNumber){
330 fLastRunNumber = fESD->GetRunNumber();
331 Prepare(fGeoFilename.Data(),fDefaultOCDB.Data(),fMisAlignOCDB.Data());
334 Int_t nTracks = Int_t(fESD->GetNumberOfMuonTracks());
335 // if (!event%100) cout << " there are " << nTracks << " tracks in event " << event << endl;
336 for (Int_t iTrack = 0; iTrack < nTracks; iTrack++) {
337 AliESDMuonTrack* esdTrack = fESD->GetMuonTrack(iTrack);
338 if (!esdTrack->ClustersStored()) continue;
339 Double_t invBenMom = esdTrack->GetInverseBendingMomentum();
340 // fInvBenMom->Fill(invBenMom);
341 // fBenMom->Fill(1./invBenMom);
342 if (TMath::Abs(invBenMom)<=1.04) {
344 AliMUONESDInterface::ESDToMUON(*esdTrack, track);
345 fAlign->ProcessTrack(&track);
346 fAlign->LocalFit(fTrackOk++,trackParams,0);
349 cout << "Processed " << fTrackTot << " Tracks." << endl;
352 // Post final data. Write histo list to a file with option "RECREATE"
356 //________________________________________________________________________
357 void AliMUONAlignmentTask::FinishTaskOutput()
359 /// Called once per task on the client machine at the end of the analysis.
361 cout << "Processed " << fTrackTot << " Tracks." << endl;
362 // Perform global fit
363 fAlign->GlobalFit(fParameters,fErrors,fPulls);
365 cout << "Done with GlobalFit " << endl;
367 // // Update pointers reading them from the output slot
368 // fList = (TList*)GetOutputData(0);
369 // fMSDEx = (TGraphErrors*)fList->At(0);
370 // fMSDEy = (TGraphErrors*)fList->At(1);
371 // fMSDEz = (TGraphErrors*)fList->At(2);
372 // fMSDEp = (TGraphErrors*)fList->At(3);
375 Double_t deId[156] = {0};
376 Double_t msdEx[156] = {0};
377 Double_t msdEy[156] = {0};
378 Double_t msdEz[156] = {0};
379 Double_t msdEp[156] = {0};
380 Double_t deIdErr[156] = {0};
381 Double_t msdExErr[156] = {0};
382 Double_t msdEyErr[156] = {0};
383 Double_t msdEzErr[156] = {0};
384 Double_t msdEpErr[156] = {0};
385 Int_t lNDetElem = 4*2+4*2+18*2+26*2+26*2;
386 Int_t lNDetElemCh[10] = {4,4,4,4,18,18,26,26,26,26};
387 // Int_t lSNDetElemCh[10] = {4,8,12,16,34,52,78,104,130,156};
388 Int_t idOffset = 0; // 400
389 Int_t lSDetElemCh = 0;
390 for(Int_t iDE=0; iDE<lNDetElem; iDE++){
392 deId[iDE] = idOffset+100;
395 for(Int_t iCh=0; iCh<9; iCh++){
396 lSDetElemCh += lNDetElemCh[iCh];
397 if (iDE>=lSDetElemCh) {
399 deId[iDE] -= lNDetElemCh[iCh];
402 msdEx[iDE]=fParameters[3*iDE+0];
403 msdEy[iDE]=fParameters[3*iDE+1];
404 msdEz[iDE]=fParameters[3*iDE+3];
405 msdEp[iDE]=fParameters[3*iDE+2];
406 msdExErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+0);
407 msdEyErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+1);
408 msdEzErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+3);
409 msdEpErr[iDE]=(Double_t)fAlign->GetParError(3*iDE+2);
410 fMSDEx->SetPoint(iDE,deId[iDE],fParameters[3*iDE+0]);
411 fMSDEx->SetPoint(iDE,deIdErr[iDE],(Double_t)fAlign->GetParError(3*iDE+0));
412 fMSDEy->SetPoint(iDE,deId[iDE],fParameters[3*iDE+1]);
413 fMSDEy->SetPoint(iDE,deIdErr[iDE],(Double_t)fAlign->GetParError(3*iDE+1));
414 fMSDEp->SetPoint(iDE,deId[iDE],fParameters[3*iDE+2]);
415 fMSDEp->SetPoint(iDE,deIdErr[iDE],(Double_t)fAlign->GetParError(3*iDE+2));
416 fMSDEz->SetPoint(iDE,deId[iDE],fParameters[3*iDE+3]);
417 fMSDEz->SetPoint(iDE,deIdErr[iDE],(Double_t)fAlign->GetParError(3*iDE+3));
420 // Post final data. Write histo list to a file with option "RECREATE"
424 AliMUONGeometryTransformer *newTransform = fAlign->ReAlign(fTransform,fParameters,true);
425 newTransform->WriteTransformations("transform2ReAlign.dat");
427 // Generate realigned data in local cdb
428 const TClonesArray* array = newTransform->GetMisAlignmentData();
430 // 100 mum residual resolution for chamber misalignments?
431 fAlign->SetAlignmentResolution(array,-1,0.01,0.01,0.004,0.003);
434 AliCDBManager* cdbManager = AliCDBManager::Instance();
435 cdbManager->SetDefaultStorage(fDefaultOCDB.Data());
436 cdbManager->SetSpecificStorage("MUON/Align/Data",fMisAlignOCDB.Data());
438 AliCDBMetaData* cdbData = new AliCDBMetaData();
439 cdbData->SetResponsible("Dimuon Offline project");
440 cdbData->SetComment("MUON alignment objects with residual misalignment");
441 AliCDBId id("MUON/Align/Data", 0, AliCDBRunRange::Infinity());
442 cdbManager->Put(const_cast<TClonesArray*>(array), id, cdbData);
446 //________________________________________________________________________
447 void AliMUONAlignmentTask::Terminate(const Option_t*)
449 /// Called once per task on the client machine at the end of the analysis.
453 //-----------------------------------------------------------------------
454 void AliMUONAlignmentTask::Prepare(const char* geoFilename, const char* defaultOCDB, const char* misAlignOCDB)
456 /// Set the geometry, the magnetic field, the mapping and the reconstruction parameters
459 AliCDBManager* man = AliCDBManager::Instance();
460 man->SetDefaultStorage(defaultOCDB);
461 man->SetSpecificStorage("MUON/Align/Data",misAlignOCDB);
463 man->SetRun(fLastRunNumber);
464 if ( ! AliMpCDB::LoadDDLStore() ) {
465 Error("MUONRefit","Could not access mapping from OCDB !");
469 // Import TGeo geometry (needed by AliMUONTrackExtrap::ExtrapToVertex)
471 AliGeomManager::LoadGeometry(geoFilename);
473 Error("AliMUONReAlignTask", "getting geometry from file %s failed", "generated/galice.root");
479 if (!TGeoGlobalMagField::Instance()->GetField()) {
480 printf("Loading field map...\n");
481 AliGRPManager *grpMan = new AliGRPManager();
482 grpMan->ReadGRPEntry();
483 grpMan->SetMagField();
486 // set the magnetic field for track extrapolations
487 AliMUONTrackExtrap::SetField();