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 **************************************************************************/
20 /// Class containing MUON data: hits, digits, rawclusters, globaltrigger, localtrigger, etc ..
21 /// The classe makes the lik between the MUON data lists and the event trees from loaders
23 /// Gines Martinez, Subatech, September 2003
26 #include "AliMUONData.h"
27 #include "AliMUONDataIterator.h"
28 #include "AliMUONConstants.h"
29 #include "AliMUONHit.h"
30 #include "AliMUONDigit.h"
31 #include "AliMUONGlobalTrigger.h"
32 #include "AliMUONLocalTrigger.h"
33 #include "AliMUONRegionalTrigger.h"
34 #include "AliMUONTriggerCrateStore.h"
35 #include "AliMUONTriggerCircuit.h"
36 #include "AliMUONGeometryTransformer.h"
37 #include "AliMUONRawCluster.h"
39 // This is from rec, classes in base should not depend on rec !!!
40 #include "AliMUONTrack.h"
41 #include "AliMUONTriggerTrack.h"
43 #include "AliRunLoader.h"
48 #include <TParticle.h>
50 #include <Riostream.h>
57 //_____________________________________________________________________________
58 AliMUONData::AliMUONData():
68 fRegionalTrigger(0x0),
70 fRecTriggerTracks(0x0),
79 fNrectriggertracks(0),
83 // Default constructor
85 //_____________________________________________________________________________
86 AliMUONData::AliMUONData(AliLoader * loader, const char* name, const char* title):
96 fRegionalTrigger(0x0),
98 fRecTriggerTracks(0x0),
105 fNregionaltrigger(0),
107 fNrectriggertracks(0),
111 /// Standard constructor
114 //_____________________________________________________________________________
115 AliMUONData::AliMUONData(const char* galiceFile):
116 TNamed("MUON", "MUON"),
125 fRegionalTrigger(0x0),
127 fRecTriggerTracks(0x0),
134 fNregionaltrigger(0),
136 fNrectriggertracks(0),
140 /// Constructor for loading data from gAlice file
142 fRunLoader = AliRunLoader::Open(galiceFile, "MUONFolder", "READ");
144 AliError(Form("Error opening %s file \n", galiceFile));
148 fLoader = fRunLoader->GetLoader("MUONLoader");
150 AliError(Form("Could get MUONLoader"));
155 //_____________________________________________________________________________
156 AliMUONData::~AliMUONData()
158 /// Destructor for AliMUONData
173 fRawClusters->Delete();
177 fGlobalTrigger->Delete();
178 delete fGlobalTrigger;
180 if (fRegionalTrigger){
181 fRegionalTrigger->Delete();
182 delete fRegionalTrigger;
185 fLocalTrigger->Delete();
186 delete fLocalTrigger;
189 fRecTracks->Delete();
192 if (fRecTriggerTracks){
193 fRecTriggerTracks->Delete();
194 delete fRecTriggerTracks;
198 fRunLoader->UnloadAll();
202 //____________________________________________________________________________
203 void AliMUONData::AddHit(Int_t fIshunt, Int_t track, Int_t detElemId,
204 Int_t idpart, Float_t X, Float_t Y, Float_t Z,
205 Float_t tof, Float_t momentum, Float_t theta,
206 Float_t phi, Float_t length, Float_t destep,
207 Float_t Xref,Float_t Yref,Float_t Zref)
209 // Add new hit to the hit list
211 TClonesArray &lhits = *fHits;
212 new(lhits[fNhits++]) AliMUONHit(fIshunt, track, detElemId,
214 tof, momentum, theta,
218 //_____________________________________________________________________________
219 void AliMUONData::AddDigit(Int_t id, const AliMUONDigit& digit)
221 /// Add a MUON digit to the list of Digits of the detection plane id
223 TClonesArray &ldigits = * Digits(id) ;
224 new(ldigits[fNdigits[id]++]) AliMUONDigit(digit);
226 //_____________________________________________________________________________
227 void AliMUONData::AddSDigit(Int_t id, const AliMUONDigit& Sdigit)
229 /// Add a MUON Sdigit to the list of SDigits of the detection plane id
231 TClonesArray &lSdigits = * SDigits(id) ;
232 new(lSdigits[fNSdigits[id]++]) AliMUONDigit(Sdigit);
235 //_____________________________________________________________________________
236 void AliMUONData::AddGlobalTrigger(const AliMUONGlobalTrigger& trigger )
238 /// Add a MUON Global Trigger to the list (only one GlobalTrigger per event !);
240 TClonesArray &globalTrigger = *fGlobalTrigger;
241 new(globalTrigger[fNglobaltrigger++]) AliMUONGlobalTrigger(trigger);
244 //____________________________________________________________________________
245 void AliMUONData::AddRegionalTrigger(const AliMUONRegionalTrigger& trigger)
247 /// add a MUON regional Trigger to the list
248 TClonesArray ®ionalTrigger = *fRegionalTrigger;
249 new(regionalTrigger[fNregionaltrigger++]) AliMUONRegionalTrigger(trigger);
251 //____________________________________________________________________________
252 void AliMUONData::AddLocalTrigger(const AliMUONLocalTrigger& trigger)
254 /// add a MUON Local Trigger to the list
256 TClonesArray &localTrigger = *fLocalTrigger;
257 new(localTrigger[fNlocaltrigger++]) AliMUONLocalTrigger(trigger);
260 //_____________________________________________________________________________
261 void AliMUONData::AddRawCluster(Int_t id, const AliMUONRawCluster& c)
263 /// Add a MUON rawcluster to the list in the detection plane id
265 TClonesArray &lrawcl = *((TClonesArray*) fRawClusters->At(id));
266 new(lrawcl[fNrawclusters[id]++]) AliMUONRawCluster(c);
268 //_____________________________________________________________________________
269 void AliMUONData::AddRecTrack(const AliMUONTrack& track)
271 /// Add a MUON rectrack
273 TClonesArray &lrectracks = *fRecTracks;
274 new(lrectracks[fNrectracks++]) AliMUONTrack(track);
276 //_____________________________________________________________________________
277 void AliMUONData::AddRecTriggerTrack(const AliMUONTriggerTrack& triggertrack)
279 /// Add a MUON triggerrectrack
281 TClonesArray &lrectriggertracks = *fRecTriggerTracks;
282 new(lrectriggertracks[fNrectriggertracks++]) AliMUONTriggerTrack(triggertrack);
284 //____________________________________________________________________________
285 TClonesArray* AliMUONData::Digits(Int_t DetectionPlane) const
287 /// Getting List of Digits
290 return ( (TClonesArray*) fDigits->At(DetectionPlane) );
294 //____________________________________________________________________________
295 TClonesArray* AliMUONData::SDigits(Int_t DetectionPlane) const
297 /// Getting List of SDigits
300 return ( (TClonesArray*) fSDigits->At(DetectionPlane) );
304 //____________________________________________________________________________
305 Bool_t AliMUONData::IsRawClusterBranchesInTree()
307 /// Checking if there are RawCluster Branches In TreeR
310 AliError("No treeR in memory");
315 sprintf(branchname,"%sRawClusters1",GetName());
316 TBranch * branch = 0x0;
317 branch = TreeR()->GetBranch(branchname);
318 if (branch) return kTRUE;
322 //____________________________________________________________________________
323 Bool_t AliMUONData::IsDigitsBranchesInTree()
325 /// Checking if there are RawCluster Branches In TreeR
328 AliError("No treeD in memory");
333 sprintf(branchname,"%sDigits1",GetName());
334 TBranch * branch = 0x0;
335 branch = TreeD()->GetBranch(branchname);
336 if (branch) return kTRUE;
340 //____________________________________________________________________________
341 Bool_t AliMUONData::IsTriggerBranchesInTree()
343 /// Checking if there are Trigger Branches In TreeR
345 AliError("No treeR in memory");
350 sprintf(branchname,"%sLocalTrigger",GetName());
351 TBranch * branch = 0x0;
352 branch = TreeR()->GetBranch(branchname);
353 if (branch) return kTRUE;
357 //____________________________________________________________________________
358 Bool_t AliMUONData::IsTriggerBranchesInTreeD()
360 /// Checking if there are Trigger Branches In TreeR
362 AliError("No treeD in memory");
367 sprintf(branchname,"%sLocalTrigger",GetName());
368 TBranch * branch = 0x0;
369 branch = TreeD()->GetBranch(branchname);
370 if (branch) return kTRUE;
375 //____________________________________________________________________________
376 Bool_t AliMUONData::IsTrackBranchesInTree()
378 /// Checking if there are Track Branches In TreeT
380 AliError("No treeT in memory");
385 sprintf(branchname,"%sTrack",GetName());
386 TBranch * branch = 0x0;
387 branch = TreeT()->GetBranch(branchname);
388 if (branch) return kTRUE;
392 //____________________________________________________________________________
393 Bool_t AliMUONData::IsTriggerTrackBranchesInTree()
395 /// Checking if there are TriggerTrack Branches In TreeT
397 AliError("No treeT in memory");
402 sprintf(branchname,"%sTriggerTrack",GetName());
403 TBranch * branch = 0x0;
404 branch = TreeT()->GetBranch(branchname);
405 if (branch) return kTRUE;
409 //____________________________________________________________________________
410 void AliMUONData::Fill(Option_t* option)
412 /// Method to fill the trees
413 const char *cH = strstr(option,"H");
414 const char *cD = strstr(option,"D"); // Digits branches in TreeD
415 const char *cS = strstr(option,"S"); // SDigits branches in TreeS
416 const char *cRC = strstr(option,"RC"); // RawCluster branches in TreeR
417 const char *cGLT = strstr(option,"GLT"); // Global and Local Trigger branches in TreeD
418 const char *cTC = strstr(option,"TC"); // global and local Trigger branches Copy in TreeR
419 const char *cRT = strstr(option,"RT"); // Reconstructed Track in TreeT
420 const char *cRL = strstr(option,"RL"); // Reconstructed Trigger Track in TreeT
422 //const char *cRP = strstr(option,"RP"); // Reconstructed Particle in TreeP
425 TBranch * branch = 0x0;
435 if ( TreeD() && cD && cGLT )
437 // Writing digits and (global+local) trigger at once.
444 if ( IsTriggerBranchesInTreeD() )
446 for (int i=0; i<AliMUONConstants::NCh(); i++)
448 sprintf(branchname,"%sDigits%d",GetName(),i+1);
449 branch = TreeD()->GetBranch(branchname);
459 if ( TreeD() && cGLT )
461 if ( IsDigitsBranchesInTree() )
463 sprintf(branchname,"%sLocalTrigger",GetName());
464 branch = TreeD()->GetBranch(branchname);
466 sprintf(branchname,"%sRegionalTrigger",GetName());
467 branch = TreeD()->GetBranch(branchname);
469 sprintf(branchname,"%sGlobalTrigger",GetName());
470 branch = TreeD()->GetBranch(branchname);
479 } // end of TreeD() handling.
489 if ( TreeR() && cRC && cTC )
495 if ( TreeR() && cRC )
497 if ( IsTriggerBranchesInTree() )
499 // Branch per branch filling
500 for (int i=0; i<AliMUONConstants::NTrackingCh(); i++)
502 sprintf(branchname,"%sRawClusters%d",GetName(),i+1);
503 branch = TreeR()->GetBranch(branchname);
515 if (IsRawClusterBranchesInTree())
517 // Branch per branch filling
518 sprintf(branchname,"%sLocalTrigger",GetName());
519 branch = TreeR()->GetBranch(branchname);
521 sprintf(branchname,"%sRegionalTrigger",GetName());
522 branch = TreeR()->GetBranch(branchname);
524 sprintf(branchname,"%sGlobalTrigger",GetName());
525 branch = TreeR()->GetBranch(branchname);
537 if ( TreeT() && cRT && cRL )
543 if ( TreeT() && cRT )
545 if (IsTriggerTrackBranchesInTree())
547 sprintf(branchname,"%sTrack",GetName());
548 branch = TreeT()->GetBranch(branchname);
557 if ( TreeT() && cRL )
559 if (IsTrackBranchesInTree())
561 sprintf(branchname,"%sTriggerTrack",GetName());
562 branch = TreeT()->GetBranch(branchname);
573 //_____________________________________________________________________________
574 void AliMUONData::MakeBranch(Option_t* option)
576 /// Create Tree branches for the MUON.
578 const Int_t kBufferSize = 4000;
581 //Setting Data Container
582 SetDataContainer(option);
584 const char *cH = strstr(option,"H");
585 const char *cD = strstr(option,"D"); // Digits branches in TreeD
586 const char *cS = strstr(option,"S"); // Digits branches in TreeS
587 const char *cRC = strstr(option,"RC"); // RawCluster branches in TreeR
588 const char *cGLT = strstr(option,"GLT"); // Global and Local Trigger branches in TreeD
589 const char *cTC = strstr(option,"TC"); // global and local Trigger branches Copy in TreeR
590 const char *cRT = strstr(option,"RT"); // Reconstructed Track in TreeT
591 const char *cRL = strstr(option,"RL"); // Reconstructed Trigger Track in TreeT
592 //const char *cRP = strstr(option,"RP"); // Reconstructed Particle in TreeP
594 TBranch * branch = 0x0;
596 // Creating Branches for Hits
598 sprintf(branchname,"%sHits",GetName());
599 branch = TreeH()->GetBranch(branchname);
601 AliInfo(Form("MakeBranch","Branch %s is already in tree.",branchname));
604 branch = TreeH()->Branch(branchname,&fHits,kBufferSize);
605 //Info("MakeBranch","Making Branch %s for hits \n",branchname);
608 //Creating Branches for Digits
617 // one branch for digits per chamber
618 for (Int_t iDetectionPlane=0; iDetectionPlane<AliMUONConstants::NCh() ;iDetectionPlane++)
620 sprintf(branchname,"%sDigits%d",GetName(),iDetectionPlane+1);
621 branch = treeD->GetBranch(branchname);
624 AliInfo(Form("Branch %s is already in tree.",branchname));
627 TClonesArray * digits = Digits(iDetectionPlane);
628 branch = treeD->Branch(branchname, &digits, kBufferSize,1);
635 // one branch for global trigger
637 sprintf(branchname,"%sGlobalTrigger",GetName());
638 branch = treeD->GetBranch(branchname);
641 AliInfo(Form("Branch GlobalTrigger is already in treeD."));
644 branch = treeD->Branch(branchname, &fGlobalTrigger, kBufferSize);
647 // one branch for regional trigger
649 sprintf(branchname,"%sRegionalTrigger",GetName());
651 branch = treeD->GetBranch(branchname);
654 AliInfo(Form("Branch RegionalTrigger is already in treeD."));
657 branch = treeD->Branch(branchname, &fRegionalTrigger, kBufferSize);
661 // one branch for local trigger
663 sprintf(branchname,"%sLocalTrigger",GetName());
665 branch = treeD->GetBranch(branchname);
668 AliInfo(Form("Branch LocalTrigger is already in treeD."));
671 branch = treeD->Branch(branchname, &fLocalTrigger, kBufferSize);
674 //Creating Branches for SDigits
675 if (TreeS() && cS ) {
676 // one branch for Sdigits per chamber
677 for (Int_t iDetectionPlane=0; iDetectionPlane<AliMUONConstants::NCh() ;iDetectionPlane++) {
678 sprintf(branchname,"%sSDigits%d",GetName(),iDetectionPlane+1);
680 branch = TreeS()->GetBranch(branchname);
682 AliInfo(Form("Branch %s is already in tree.",branchname));
685 TClonesArray * sdigits = SDigits(iDetectionPlane);
686 branch = TreeS()->Branch(branchname, &sdigits, kBufferSize,1);
687 //Info("MakeBranch","Making Branch %s for sdigits in detection plane %d\n",branchname,iDetectionPlane+1);
691 if (TreeR() && cRC ) {
692 // one branch for raw clusters per tracking detection plane
695 for (i=0; i<AliMUONConstants::NTrackingCh() ;i++) {
696 sprintf(branchname,"%sRawClusters%d",GetName(),i+1);
698 branch = TreeR()->GetBranch(branchname);
700 AliInfo(Form("Branch %s is already in tree.",branchname));
703 branch = TreeR()->Branch(branchname, &((*fRawClusters)[i]),kBufferSize);
704 //Info("MakeBranch","Making Branch %s for rawcluster in detection plane %d\n",branchname,i+1);
708 if (TreeR() && cTC ) {
710 // one branch for global trigger
712 sprintf(branchname,"%sGlobalTrigger",GetName());
714 branch = TreeR()->GetBranch(branchname);
716 AliInfo(Form("Branch GlobalTrigger is already in treeR."));
719 branch = TreeR()->Branch(branchname, &fGlobalTrigger, kBufferSize);
720 //Info("MakeBranch", "Making Branch %s for Global Trigger\n",branchname);
723 // one branch for regional trigger
725 sprintf(branchname,"%sRegionalTrigger",GetName());
727 branch = TreeR()->GetBranch(branchname);
729 AliInfo(Form("Branch RegionalTrigger is already in treeR."));
732 branch = TreeR()->Branch(branchname, &fRegionalTrigger, kBufferSize);
735 // one branch for local trigger
737 sprintf(branchname,"%sLocalTrigger",GetName());
739 branch = TreeR()->GetBranch(branchname);
741 AliInfo(Form("Branch LocalTrigger is already in treeR."));
744 branch = TreeR()->Branch(branchname, &fLocalTrigger, kBufferSize);
745 //Info("MakeBranch", "Making Branch %s for Global Trigger\n",branchname);
748 if (TreeT() && cRT ) {
749 sprintf(branchname,"%sTrack",GetName());
750 branch = TreeT()->GetBranch(branchname);
752 AliInfo(Form("Branch %s is already in tree.",GetName()));
755 branch = TreeT()->Branch(branchname,&fRecTracks,kBufferSize);
756 //Info("MakeBranch","Making Branch %s for tracks \n",branchname);
759 if (TreeT() && cRL ) {
760 sprintf(branchname,"%sTriggerTrack",GetName());
761 branch = TreeT()->GetBranch(branchname);
763 AliInfo(Form("Branch %s is already in tree.",GetName()));
766 branch = TreeT()->Branch(branchname,&fRecTriggerTracks,kBufferSize);
767 //Info("MakeBranch","Making Branch %s for trigger tracks \n",branchname);
770 //____________________________________________________________________________
771 TClonesArray* AliMUONData::RawClusters(Int_t DetectionPlane)
773 /// Getting Raw Clusters
776 return ( (TClonesArray*) fRawClusters->At(DetectionPlane) );
781 //____________________________________________________________________________
783 AliMUONData::LocalTrigger() const
785 /// Getting local trigger
787 return fLocalTrigger;
790 //____________________________________________________________________________
792 AliMUONData::RegionalTrigger() const
794 /// Getting regional trigger
796 return fRegionalTrigger;
799 //____________________________________________________________________________
801 AliMUONData::GetNtracks() const
803 /// Get number of entries in hits three
806 if (fLoader && fLoader->TreeH())
807 ntrk = (Int_t) fLoader->TreeH()->GetEntries();
811 //____________________________________________________________________________
813 AliMUONData::GetDigits() const
815 /// Load the digits from TreeD for the current event.
817 Int_t event = fLoader->GetRunLoader()->GetEventNumber();
818 if ( fCurrentEvent != event )
820 if (fLoader->TreeD()) {
821 fLoader->TreeD()->GetEvent(0);
822 fCurrentEvent = event;
827 //____________________________________________________________________________
829 AliMUONData::GlobalTrigger() const
831 /// Return the global trigger
833 return fGlobalTrigger;
836 //____________________________________________________________________________
837 void AliMUONData::ResetDigits()
839 /// Reset number of digits and the digits array for this detector
841 if (fDigits == 0x0) return;
842 for ( int i=0;i<AliMUONConstants::NCh();i++ ) {
843 if ((*fDigits)[i]) ((TClonesArray*)fDigits->At(i))->Clear("C");
844 if (fNdigits) fNdigits[i]=0;
847 //____________________________________________________________________________
848 void AliMUONData::ResetSDigits()
850 /// Reset number of Sdigits and the Sdigits array for this detector
852 if (fSDigits == 0x0) return;
853 for ( int i=0;i<AliMUONConstants::NCh();i++ ) {
854 if ((*fSDigits)[i]) ((TClonesArray*)fSDigits->At(i))->Clear();
855 if (fNSdigits) fNSdigits[i]=0;
858 //______________________________________________________________________________
859 void AliMUONData::ResetHits()
861 /// Reset number of clusters and the cluster array for this detector
864 if (fHits) fHits->Clear();
866 //_______________________________________________________________________________
867 void AliMUONData::ResetRawClusters()
869 /// Reset number of raw clusters and the raw clust array for this detector
871 for ( int i=0;i<AliMUONConstants::NTrackingCh();i++ ) {
872 if ((*fRawClusters)[i]) ((TClonesArray*)fRawClusters->At(i))->Clear();
873 if (fNrawclusters) fNrawclusters[i]=0;
876 //_______________________________________________________________________________
877 void AliMUONData::ResetTrigger()
879 /// Reset Local and Global Trigger
882 if (fGlobalTrigger) fGlobalTrigger->Clear();
883 fNregionaltrigger = 0;
884 if (fRegionalTrigger) fRegionalTrigger->Clear();
886 if (fLocalTrigger) fLocalTrigger->Clear();
889 //____________________________________________________________________________
890 void AliMUONData::ResetRecTracks()
892 /// Reset tracks information
895 if (fRecTracks) fRecTracks->Delete(); // necessary to delete in case of memory allocation
897 //____________________________________________________________________________
898 void AliMUONData::ResetRecTriggerTracks()
900 /// Reset tracks information
902 fNrectriggertracks = 0;
903 if (fRecTriggerTracks) fRecTriggerTracks->Delete(); // necessary to delete in case of memory allocation
905 //____________________________________________________________________________
906 void AliMUONData::SetDataContainer(Option_t* option)
908 /// Setting data containers of muon data
909 const char *cH = strstr(option,"H");
910 const char *cD = strstr(option,"D"); // Digits
911 const char *cS = strstr(option,"S"); // SDigits
912 const char *cRC = strstr(option,"RC"); // RawCluster
913 const char *cGLT = strstr(option,"GLT"); // Global and Local Trigger
914 const char *cTC = strstr(option,"TC"); // global and local Trigger
915 const char *cRT = strstr(option,"RT"); // Reconstructed Tracks
916 const char *cRL = strstr(option,"RL"); // Reconstructed Trigger Tracks
917 //const char *cRP = strstr(option,"RP"); // Reconstructed Particles
918 AliDebug(1,Form("option=%s",option));
920 // Clones array for hits
923 fHits = new TClonesArray("AliMUONHit",1000);
929 // ObjArray of ClonesArrays for Digits
931 if (fDigits == 0x0 ) {
932 fDigits = new TObjArray(AliMUONConstants::NCh());
933 fNdigits= new Int_t[AliMUONConstants::NCh()];
934 for (Int_t i=0; i<AliMUONConstants::NCh() ;i++) {
935 TClonesArray * tca = new TClonesArray("AliMUONDigit",10000);
937 fDigits->AddAt(tca,i);
942 AliDebug(1,Form("fDigits already there = %p",fSDigits));
948 // ClonesArrays for Trigger
950 if (fLocalTrigger == 0x0) {
951 fLocalTrigger = new TClonesArray("AliMUONLocalTrigger",234);
953 if (fRegionalTrigger == 0x0) {
954 fRegionalTrigger = new TClonesArray("AliMUONRegionalTrigger",16);
956 if (fGlobalTrigger== 0x0) {
957 fGlobalTrigger = new TClonesArray("AliMUONGlobalTrigger",1);
963 // Container for Sdigits
965 if (fSDigits == 0x0) {
966 AliDebug(1,"Creating fSDigits TObjArray");
967 fSDigits = new TObjArray(AliMUONConstants::NCh());
968 fNSdigits= new Int_t[AliMUONConstants::NCh()];
969 for (Int_t i=0; i<AliMUONConstants::NCh() ;i++) {
970 TClonesArray* a = new TClonesArray("AliMUONDigit",10000);
972 fSDigits->AddAt(a,i);
973 AliDebug(1,Form("fSDigits[%d]=%p",i,a));
978 AliDebug(1,Form("fSDigits already there = %p",fSDigits));
984 // Containers for rawclusters, globaltrigger and local trigger tree
986 if (fRawClusters == 0x0) {
987 fRawClusters = new TObjArray(AliMUONConstants::NTrackingCh());
988 fNrawclusters= new Int_t[AliMUONConstants::NTrackingCh()];
989 for (Int_t i=0; i<AliMUONConstants::NTrackingCh();i++) {
990 TClonesArray* tca = new TClonesArray("AliMUONRawCluster",10000);
992 fRawClusters->AddAt(tca,i);
996 // ResetRawClusters();
997 // It breaks the correct functioning of the combined reconstruction (AZ)
1001 if (fLocalTrigger == 0x0) {
1002 fLocalTrigger = new TClonesArray("AliMUONLocalTrigger",234);
1004 if (fRegionalTrigger == 0x0) {
1005 fRegionalTrigger = new TClonesArray("AliMUONRegionalTrigger",16);
1007 if (fGlobalTrigger== 0x0) {
1008 fGlobalTrigger = new TClonesArray("AliMUONGlobalTrigger",1);
1011 // This is not necessary here since trigger info ins copied from digits info on flight to RecPoint output
1015 // Containers for rectracks and rectrigger tracks
1017 if (fRecTracks == 0x0) {
1018 fRecTracks = new TClonesArray("AliMUONTrack",100);
1023 if (fRecTriggerTracks == 0x0 && cRL) {
1024 fRecTriggerTracks = new TClonesArray("AliMUONTriggerTrack",100);
1026 ResetRecTriggerTracks();
1030 //____________________________________________________________________________
1031 void AliMUONData::SetTreeAddress(Option_t* option)
1033 // Setting Data containers
1034 SetDataContainer(option);
1036 /// Setting Addresses to the events trees
1038 const char *cH = strstr(option,"H");
1039 const char *cD = strstr(option,"D"); // Digits branches in TreeD
1040 const char *cS = strstr(option,"S"); // SDigits branches in TreeS
1041 const char *cRC = strstr(option,"RC"); // RawCluster branches in TreeR
1042 const char *cGLT = strstr(option,"GLT"); // Global and Local Trigger branches in TreeD
1043 const char *cTC = strstr(option,"TC"); // global and local Trigger branches Copy in TreeR
1044 const char *cRT = strstr(option,"RT"); // Reconstructed Track in TreeT
1045 const char *cRL = strstr(option,"RL"); // Reconstructed Trigger Track in TreeT
1046 //const char *cRP = strstr(option,"RP"); // Reconstructed Particle in TreeP
1048 // Set branch address for the Hits, Digits, RawClusters, GlobalTrigger and LocalTrigger Tree.
1049 char branchname[30];
1050 TBranch * branch = 0x0;
1052 AliDebug(1,Form("option=%s",option));
1054 // Branch address for hit tree
1055 if (TreeH() && fHits && cH) {
1056 sprintf(branchname,"%sHits",GetName());
1057 branch = TreeH()->GetBranch(branchname);
1059 // Info("SetTreeAddress","(%s) Setting for Hits",GetName());
1060 branch->SetAddress(&fHits);
1062 else { //can be invoked before branch creation
1063 //AliWarning(Form("(%s) Failed for Hits. Can not find branch in tree.",GetName()));
1068 // Branch address for digit tree
1069 if (TreeD() && fDigits && cD) {
1070 for (int i=0; i<AliMUONConstants::NCh(); i++) {
1071 sprintf(branchname,"%sDigits%d",GetName(),i+1);
1073 branch = TreeD()->GetBranch(branchname);
1074 TClonesArray * digits = Digits(i);
1076 branch->SetAddress( &digits );
1078 else AliWarning(Form("(%s) Failed for Digits Detection plane %d. Can not find branch in tree.",GetName(),i));
1082 if ( TreeD() && fLocalTrigger && cGLT) {
1083 sprintf(branchname,"%sLocalTrigger",GetName());
1084 branch = TreeD()->GetBranch(branchname);
1085 if (branch) branch->SetAddress(&fLocalTrigger);
1086 else AliWarning(Form("(%s) Failed for LocalTrigger. Can not find branch in treeD.",GetName()));
1088 if ( TreeD() && fRegionalTrigger && cGLT) {
1089 sprintf(branchname,"%sRegionalTrigger",GetName());
1090 branch = TreeD()->GetBranch(branchname);
1091 if (branch) branch->SetAddress(&fRegionalTrigger);
1092 else AliWarning(Form("(%s) Failed for RegionalTrigger. Can not find branch in treeD.",GetName()));
1094 if ( TreeD() && fGlobalTrigger && cGLT) {
1095 sprintf(branchname,"%sGlobalTrigger",GetName());
1096 branch = TreeD()->GetBranch(branchname);
1097 if (branch) branch->SetAddress(&fGlobalTrigger);
1098 else AliWarning(Form("(%s) Failed for GlobalTrigger. Can not find branch in treeD.",GetName()));
1102 // Branch address for Sdigit tree
1103 if (TreeS() && fSDigits && cS) {
1104 AliDebug(1,"Setting branch addresses");
1105 for (int i=0; i<AliMUONConstants::NCh(); i++) {
1106 sprintf(branchname,"%sSDigits%d",GetName(),i+1);
1108 AliDebug(1,Form("TreeS=%p for ich=%d branchname=%s",
1109 TreeS(),i,branchname));
1110 branch = TreeS()->GetBranch(branchname);
1111 TClonesArray * sdigits = SDigits(i);
1112 if (branch) branch->SetAddress( &sdigits );
1113 else AliWarning(Form("(%s) Failed for SDigits Detection plane %d. Can not find branch in tree.",GetName(),i));
1119 // Branch address for rawclusters, globaltrigger and local trigger tree
1120 if ( TreeR() && fRawClusters && cRC && !strstr(cRC,"RCC")) {
1121 for (int i=0; i<AliMUONConstants::NTrackingCh(); i++) {
1122 sprintf(branchname,"%sRawClusters%d",GetName(),i+1);
1124 branch = TreeR()->GetBranch(branchname);
1125 if (branch) branch->SetAddress( &((*fRawClusters)[i]) );
1126 else AliWarning(Form("(%s) Failed for RawClusters Detection plane %d. Can not find branch in tree.",GetName(),i));
1130 if ( TreeR() && fLocalTrigger && cTC) {
1131 sprintf(branchname,"%sLocalTrigger",GetName());
1132 branch = TreeR()->GetBranch(branchname);
1133 if (branch) branch->SetAddress(&fLocalTrigger);
1134 else AliWarning(Form("(%s) Failed for LocalTrigger. Can not find branch in treeR.",GetName()));
1136 if ( TreeR() && fRegionalTrigger && cTC) {
1137 sprintf(branchname,"%sRegionalTrigger",GetName());
1138 branch = TreeR()->GetBranch(branchname);
1139 if (branch) branch->SetAddress(&fRegionalTrigger);
1140 else AliWarning(Form("(%s) Failed for RegionalTrigger. Can not find branch in treeR.",GetName()));
1142 if ( TreeR() && fGlobalTrigger && cTC) {
1143 sprintf(branchname,"%sGlobalTrigger",GetName());
1144 branch = TreeR()->GetBranch(branchname);
1145 if (branch) branch->SetAddress(&fGlobalTrigger);
1146 else AliWarning(Form("(%s) Failed for GlobalTrigger. Can not find branch in treeR.",GetName()));
1150 if ( TreeT() && fRecTracks && cRT ) {
1151 sprintf(branchname,"%sTrack",GetName());
1152 branch = TreeT()->GetBranch(branchname);
1153 if (branch) branch->SetAddress(&fRecTracks);
1154 else AliWarning(Form("(%s) Failed for Tracks. Can not find branch in tree.",GetName()));
1157 if ( TreeT() && fRecTriggerTracks && cRL ) {
1158 sprintf(branchname,"%sTriggerTrack",GetName());
1159 branch = TreeT()->GetBranch(branchname);
1160 if (branch) branch->SetAddress(&fRecTriggerTracks);
1161 else AliWarning(Form("(%s) Failed for Trigger Tracks. Can not find branch in tree.",GetName()));
1165 //_____________________________________________________________________________
1167 AliMUONData::Print(Option_t* opt) const
1169 /// Dump object on screen
1171 TString options(opt);
1174 if ( options.Contains("D") )
1176 for ( Int_t ich = 0; ich < AliMUONConstants::NCh(); ++ich)
1178 TClonesArray* digits = Digits(ich);
1179 Int_t ndigits = digits->GetEntriesFast();
1180 for ( Int_t id = 0; id < ndigits; ++id )
1182 AliMUONDigit* digit =
1183 static_cast<AliMUONDigit*>(digits->UncheckedAt(id));
1189 if ( options.Contains("S") )
1191 for ( Int_t ich = 0; ich < AliMUONConstants::NCh(); ++ich)
1193 TClonesArray* digits = SDigits(ich);
1194 Int_t ndigits = digits->GetEntriesFast();
1195 for ( Int_t id = 0; id < ndigits; ++id )
1197 AliMUONDigit* digit =
1198 static_cast<AliMUONDigit*>(digits->UncheckedAt(id));
1206 //_____________________________________________________________________________
1208 AliMUONData::DumpKine(Int_t event2Check)
1211 fRunLoader->LoadKinematics("READ");
1213 Int_t nevents = fRunLoader->GetNumberOfEvents();
1214 for (Int_t ievent=0; ievent<nevents; ievent++) { // Event loop
1215 if ( event2Check != 0 ) ievent=event2Check;
1217 // Getting event ievent
1218 fRunLoader->GetEvent(ievent);
1220 // Stack of particle for this event
1221 AliStack* stack = fRunLoader->Stack();
1223 Int_t nparticles = (Int_t) fRunLoader->Stack()->GetNtrack();
1224 printf(">>> Event %d, Number of particles is %d \n", ievent, nparticles);
1226 for (Int_t iparticle=0; iparticle<nparticles; iparticle++) {
1227 stack->Particle(iparticle)->Print("");
1229 if (event2Check!=0) ievent=nevents;
1231 fRunLoader->UnloadKinematics();
1235 //_____________________________________________________________________________
1237 AliMUONData::DumpHits(Int_t event2Check, Option_t* opt)
1240 fLoader->LoadHits("READ");
1243 Int_t nevents = fRunLoader->GetNumberOfEvents();
1244 for (Int_t ievent=0; ievent<nevents; ievent++) {
1245 if (event2Check!=0) ievent=event2Check;
1246 printf(">>> Event %d \n",ievent);
1248 // Getting event ievent
1249 fRunLoader->GetEvent(ievent);
1250 SetTreeAddress("H");
1253 Int_t ntracks = (Int_t) GetNtracks();
1254 for (Int_t itrack=0; itrack<ntracks; itrack++) {
1255 //Getting List of Hits of Track itrack
1258 Int_t nhits = (Int_t) Hits()->GetEntriesFast();
1259 printf(">>> Track %d, Number of hits %d \n",itrack,nhits);
1260 for (Int_t ihit=0; ihit<nhits; ihit++) {
1261 AliMUONHit* mHit = static_cast<AliMUONHit*>(Hits()->At(ihit));
1266 if (event2Check!=0) ievent=nevents;
1268 fLoader->UnloadHits();
1271 //_____________________________________________________________________________
1273 AliMUONData::DumpDigits(Int_t event2Check, Option_t* opt)
1276 fLoader->LoadDigits("READ");
1279 Int_t firstEvent = 0;
1280 Int_t lastEvent = fRunLoader->GetNumberOfEvents()-1;
1281 if ( event2Check != 0 ) {
1282 firstEvent = event2Check;
1283 lastEvent = event2Check;
1286 for ( Int_t ievent = firstEvent; ievent <= lastEvent; ++ievent ) {
1287 printf(">>> Event %d \n",ievent);
1288 fRunLoader->GetEvent(ievent);
1290 AliMUONDataIterator it(this, "digit", AliMUONDataIterator::kTrackingChambers);
1291 AliMUONDigit* digit;
1293 while ( ( digit = (AliMUONDigit*)it.Next() ) )
1298 fLoader->UnloadDigits();
1301 //_____________________________________________________________________________
1303 AliMUONData::DumpSDigits(Int_t event2Check, Option_t* opt)
1306 fLoader->LoadSDigits("READ");
1309 Int_t nevents = fRunLoader->GetNumberOfEvents();
1310 for (Int_t ievent=0; ievent<nevents; ievent++) {
1311 if (event2Check!=0) ievent=event2Check;
1312 printf(">>> Event %d \n",ievent);
1314 // Getting event ievent
1315 fRunLoader->GetEvent(ievent);
1316 SetTreeAddress("S");
1320 Int_t nchambers = AliMUONConstants::NCh(); ;
1321 for (Int_t ichamber=0; ichamber<nchambers; ichamber++) {
1322 TClonesArray* digits = SDigits(ichamber);
1325 Int_t ndigits = (Int_t)digits->GetEntriesFast();
1326 for (Int_t idigit=0; idigit<ndigits; idigit++) {
1327 AliMUONDigit* mDigit = static_cast<AliMUONDigit*>(digits->At(idigit));
1332 if (event2Check!=0) ievent=nevents;
1334 fLoader->UnloadSDigits();
1337 //_____________________________________________________________________________
1339 AliMUONData::DumpRecPoints(Int_t event2Check, Option_t* opt)
1342 fLoader->LoadRecPoints("READ");
1345 Int_t nevents = fRunLoader->GetNumberOfEvents();
1346 for (Int_t ievent=0; ievent<nevents; ievent++) {
1347 if (event2Check!=0) ievent=event2Check;
1348 printf(">>> Event %d \n",ievent);
1350 // Getting event ievent
1351 fRunLoader->GetEvent(ievent);
1352 Int_t nchambers = AliMUONConstants::NTrackingCh();
1353 SetTreeAddress("RC,TC");
1357 for (Int_t ichamber=0; ichamber<nchambers; ichamber++) {
1358 char branchname[30];
1359 sprintf(branchname,"MUONRawClusters%d",ichamber+1);
1360 //printf(">>> branchname %s\n",branchname);
1362 // Loop on rec points
1363 Int_t nrecpoints = (Int_t) RawClusters(ichamber)->GetEntriesFast();
1364 // printf(">>> Chamber %2d, Number of recpoints = %6d \n",ichamber+1, nrecpoints);
1365 for (Int_t irecpoint=0; irecpoint<nrecpoints; irecpoint++) {
1366 AliMUONRawCluster* mRecPoint = static_cast<AliMUONRawCluster*>(RawClusters(ichamber)->At(irecpoint));
1367 mRecPoint->Print(opt);
1371 if (event2Check!=0) ievent=nevents;
1373 fLoader->UnloadRecPoints();
1377 //_____________________________________________________________________________
1379 AliMUONData::DumpRecTrigger(Int_t event2Check,
1380 Int_t write, Bool_t readFromRP)
1382 /// Reads and dumps trigger objects from MUON.RecPoints.root
1384 TClonesArray * globalTrigger;
1385 TClonesArray * localTrigger;
1387 // Do NOT print out all the info if the loop runs over all events
1388 Int_t printout = (event2Check == 0 ) ? 0 : 1 ;
1390 // Book a ntuple for more detailled studies
1391 TNtuple *tupleGlo = new TNtuple("TgtupleGlo","Global Trigger Ntuple","ev:global:slpt:shpt:uplpt:uphpt:lplpt:lplpt");
1392 TNtuple *tupleLoc = new TNtuple("TgtupleLoc","Local Trigger Ntuple","ev:LoCircuit:LoStripX:LoDev:StripY:LoLpt:LoHpt:y11:y21:x11");
1395 Int_t sLowpt=0,sHighpt=0;
1396 Int_t uSLowpt=0,uSHighpt=0;
1397 Int_t lSLowpt=0,lSHighpt=0;
1399 AliMUONTriggerCrateStore* crateManager = new AliMUONTriggerCrateStore();
1400 crateManager->ReadFromFile();
1402 AliMUONGeometryTransformer* transformer = new AliMUONGeometryTransformer(kFALSE);
1403 transformer->ReadGeometryData("volpath.dat", "geometry.root");
1405 TClonesArray* triggerCircuit = new TClonesArray("AliMUONTriggerCircuit", 234);
1407 for (Int_t i = 0; i < AliMUONConstants::NTriggerCircuit(); i++) {
1408 AliMUONTriggerCircuit* c = new AliMUONTriggerCircuit();
1409 c->SetTransformer(transformer);
1410 c->Init(i,*crateManager);
1411 TClonesArray& circuit = *triggerCircuit;
1412 new(circuit[circuit.GetEntriesFast()])AliMUONTriggerCircuit(*c);
1416 Char_t fileName[30];
1418 AliInfoStream() << " reading from digits \n";
1419 fLoader->LoadDigits("READ");
1420 sprintf(fileName,"TriggerCheckFromDigits.root");
1422 AliInfoStream() << " reading from RecPoints \n";
1423 fLoader->LoadRecPoints("READ");
1424 sprintf(fileName,"TriggerCheckFromRP.root");
1428 AliMUONGlobalTrigger *gloTrg(0x0);
1429 AliMUONLocalTrigger *locTrg(0x0);
1431 Int_t nevents = fRunLoader->GetNumberOfEvents();
1432 for (Int_t ievent=0; ievent<nevents; ievent++) {
1433 if (event2Check!=0) ievent=event2Check;
1434 if (ievent%100==0 || event2Check)
1435 AliInfoStream() << "Processing event " << ievent << endl;
1436 fRunLoader->GetEvent(ievent);
1439 SetTreeAddress("D,GLT");
1442 SetTreeAddress("RC,TC");
1446 globalTrigger = GlobalTrigger();
1447 localTrigger = LocalTrigger();
1449 Int_t nglobals = (Int_t) globalTrigger->GetEntriesFast(); // should be 1
1450 Int_t nlocals = (Int_t) localTrigger->GetEntriesFast(); // up to 234
1451 if (printout) printf("###################################################\n");
1452 if (printout) printf("event %d nglobal %d nlocal %d \n",ievent,nglobals,nlocals);
1454 for (Int_t iglobal=0; iglobal<nglobals; iglobal++) { // Global Trigger
1455 gloTrg = static_cast<AliMUONGlobalTrigger*>(globalTrigger->At(iglobal));
1457 sLowpt+=gloTrg->SingleLpt() ;
1458 sHighpt+=gloTrg->SingleHpt() ;
1459 uSLowpt+=gloTrg->PairUnlikeLpt();
1460 uSHighpt+=gloTrg->PairUnlikeHpt();
1461 lSLowpt+=gloTrg->PairLikeLpt();
1462 lSHighpt+=gloTrg->PairLikeHpt();
1464 if (printout) gloTrg->Print("full");
1466 } // end of loop on Global Trigger
1468 for (Int_t ilocal=0; ilocal<nlocals; ilocal++) { // Local Trigger
1469 locTrg = static_cast<AliMUONLocalTrigger*>(localTrigger->At(ilocal));
1471 if (locTrg->LoLpt()!=0) { // board is fired
1473 if (printout) locTrg->Print("full");
1475 AliMUONTriggerCircuit* circuit = (AliMUONTriggerCircuit*)triggerCircuit->At(locTrg->LoCircuit()-1);
1477 tupleLoc->Fill(ievent,locTrg->LoCircuit(),locTrg->LoStripX(),locTrg->LoDev(),locTrg->LoStripY(),locTrg->LoLpt(),locTrg->LoHpt(),circuit->GetY11Pos(locTrg->LoStripX()),circuit->GetY21Pos(locTrg->LoStripX()+locTrg->LoDev()+1),circuit->GetX11Pos(locTrg->LoStripY()));
1480 } // end of loop on Local Trigger
1483 tupleGlo->Fill(ievent,nglobals,gloTrg->SingleLpt(),gloTrg->SingleHpt(),gloTrg->PairUnlikeLpt(),gloTrg->PairUnlikeHpt(),gloTrg->PairLikeLpt(),gloTrg->PairLikeHpt());
1486 if (event2Check!=0) ievent=nevents;
1487 } // end loop on event
1489 // Print out summary if loop ran over all event
1493 printf("=============================================\n");
1494 printf("================ SUMMARY ==================\n");
1496 printf("Total number of events processed %d \n", (event2Check==0) ? nevents : 1);
1498 printf(" Global Trigger output Low pt High pt\n");
1499 printf(" number of Single :\t");
1500 printf("%i\t%i\t",sLowpt,sHighpt);
1502 printf(" number of UnlikeSign pair :\t");
1503 printf("%i\t%i\t",uSLowpt,uSHighpt);
1505 printf(" number of LikeSign pair :\t");
1506 printf("%i\t%i\t",lSLowpt,lSHighpt);
1508 printf("=============================================\n");
1513 TFile *myFile = new TFile(fileName, "RECREATE");
1519 fLoader->UnloadRecPoints();
1521 delete crateManager;
1523 delete triggerCircuit;