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 Revision 1.18 2000/06/12 15:15:46 jbarbosa
21 Revision 1.17 2000/06/09 14:58:37 jbarbosa
22 New digitisation per particle type
24 Revision 1.16 2000/04/19 12:55:43 morsch
25 Newly structured and updated version (JB, AM)
30 ////////////////////////////////////////////////
31 // Manager and hits classes for set:RICH //
32 ////////////////////////////////////////////////
40 #include <TObjArray.h>
43 #include <TParticle.h>
47 #include "AliRICHSegmentation.h"
48 #include "AliRICHHit.h"
49 #include "AliRICHCerenkov.h"
50 #include "AliRICHPadHit.h"
51 #include "AliRICHDigit.h"
52 #include "AliRICHTransientDigit.h"
53 #include "AliRICHRawCluster.h"
54 #include "AliRICHRecHit.h"
55 #include "AliRICHHitMapA1.h"
56 #include "AliRICHClusterFinder.h"
59 #include "AliPoints.h"
60 #include "AliCallf77.h"
63 // Static variables for the pad-hit iterator routines
64 static Int_t sMaxIterPad=0;
65 static Int_t sCurIterPad=0;
66 static TClonesArray *fClusters2;
67 static TClonesArray *fHits2;
72 //___________________________________________
75 // Default constructor for RICH manager class
87 //___________________________________________
88 AliRICH::AliRICH(const char *name, const char *title)
89 : AliDetector(name,title)
93 <img src="gif/alirich.gif">
97 fHits = new TClonesArray("AliRICHHit",1000 );
98 gAlice->AddHitList(fHits);
99 fPadHits = new TClonesArray("AliRICHPadHit",100000);
100 fCerenkovs = new TClonesArray("AliRICHCerenkov",1000);
101 gAlice->AddHitList(fCerenkovs);
102 //gAlice->AddHitList(fHits);
107 fNdch = new Int_t[kNCH];
109 fDchambers = new TObjArray(kNCH);
111 fRecHits = new TObjArray(kNCH);
115 for (i=0; i<kNCH ;i++) {
116 (*fDchambers)[i] = new TClonesArray("AliRICHDigit",10000);
120 fNrawch = new Int_t[kNCH];
122 fRawClusters = new TObjArray(kNCH);
123 //printf("Created fRwClusters with adress:%p",fRawClusters);
125 for (i=0; i<kNCH ;i++) {
126 (*fRawClusters)[i] = new TClonesArray("AliRICHRawCluster",10000);
130 fNrechits = new Int_t[kNCH];
132 for (i=0; i<kNCH ;i++) {
133 (*fRecHits)[i] = new TClonesArray("AliRICHRecHit",1000);
135 //printf("Created fRecHits with adress:%p",fRecHits);
138 SetMarkerColor(kRed);
141 AliRICH::AliRICH(const AliRICH& RICH)
147 //___________________________________________
151 // Destructor of RICH manager class
159 //___________________________________________
160 void AliRICH::AddHit(Int_t track, Int_t *vol, Float_t *hits)
164 // Adds a hit to the Hits list
167 TClonesArray &lhits = *fHits;
168 new(lhits[fNhits++]) AliRICHHit(fIshunt,track,vol,hits);
170 //_____________________________________________________________________________
171 void AliRICH::AddCerenkov(Int_t track, Int_t *vol, Float_t *cerenkovs)
175 // Adds a RICH cerenkov hit to the Cerenkov Hits list
178 TClonesArray &lcerenkovs = *fCerenkovs;
179 new(lcerenkovs[fNcerenkovs++]) AliRICHCerenkov(fIshunt,track,vol,cerenkovs);
180 //printf ("Done for Cerenkov %d\n\n\n\n",fNcerenkovs);
182 //___________________________________________
183 void AliRICH::AddPadHit(Int_t *clhits)
187 // Add a RICH pad hit to the list
190 TClonesArray &lPadHits = *fPadHits;
191 new(lPadHits[fNPadHits++]) AliRICHPadHit(clhits);
193 //_____________________________________________________________________________
194 void AliRICH::AddDigits(Int_t id, Int_t *tracks, Int_t *charges, Int_t *digits)
198 // Add a RICH digit to the list
201 TClonesArray &ldigits = *((TClonesArray*)(*fDchambers)[id]);
202 new(ldigits[fNdch[id]++]) AliRICHDigit(tracks,charges,digits);
205 //_____________________________________________________________________________
206 void AliRICH::AddRawCluster(Int_t id, const AliRICHRawCluster& c)
209 // Add a RICH digit to the list
212 TClonesArray &lrawcl = *((TClonesArray*)(*fRawClusters)[id]);
213 new(lrawcl[fNrawch[id]++]) AliRICHRawCluster(c);
216 //_____________________________________________________________________________
217 void AliRICH::AddRecHit(Int_t id, Float_t *rechit, Float_t *photons, Int_t *padsx, Int_t* padsy)
221 // Add a RICH reconstructed hit to the list
224 TClonesArray &lrec = *((TClonesArray*)(*fRecHits)[id]);
225 new(lrec[fNrechits[id]++]) AliRICHRecHit(id,rechit,photons,padsx,padsy);
228 //___________________________________________
229 void AliRICH::BuildGeometry()
234 // Builds a TNode geometry for event display
238 const int kColorRICH = kGreen;
240 top=gAlice->GetGeometry()->GetNode("alice");
243 new TBRIK("S_RICH","S_RICH","void",71.09999,11.5,73.15);
246 Float_t pos1[3]={0,471.8999,165.2599};
247 //Chamber(0).SetChamberTransform(pos1[0],pos1[1],pos1[2],
248 new TRotMatrix("rot993","rot993",90,0,70.69,90,19.30999,-90);
249 node = new TNode("RICH1","RICH1","S_RICH",pos1[0],pos1[1],pos1[2],"rot993");
252 node->SetLineColor(kColorRICH);
256 Float_t pos2[3]={171,470,0};
257 //Chamber(1).SetChamberTransform(pos2[0],pos2[1],pos2[2],
258 new TRotMatrix("rot994","rot994",90,-20,90,70,0,0);
259 node = new TNode("RICH2","RICH2","S_RICH",pos2[0],pos2[1],pos2[2],"rot994");
262 node->SetLineColor(kColorRICH);
265 Float_t pos3[3]={0,500,0};
266 //Chamber(2).SetChamberTransform(pos3[0],pos3[1],pos3[2],
267 new TRotMatrix("rot995","rot995",90,0,90,90,0,0);
268 node = new TNode("RICH3","RICH3","S_RICH",pos3[0],pos3[1],pos3[2],"rot995");
271 node->SetLineColor(kColorRICH);
274 Float_t pos4[3]={-171,470,0};
275 //Chamber(3).SetChamberTransform(pos4[0],pos4[1],pos4[2],
276 new TRotMatrix("rot996","rot996",90,20,90,110,0,0);
277 node = new TNode("RICH4","RICH4","S_RICH",pos4[0],pos4[1],pos4[2],"rot996");
280 node->SetLineColor(kColorRICH);
283 Float_t pos5[3]={161.3999,443.3999,-165.3};
284 //Chamber(4).SetChamberTransform(pos5[0],pos5[1],pos5[2],
285 new TRotMatrix("rot997","rot997",90,340,108.1999,70,18.2,70);
286 node = new TNode("RICH5","RICH5","S_RICH",pos5[0],pos5[1],pos5[2],"rot997");
288 node->SetLineColor(kColorRICH);
291 Float_t pos6[3]={0., 471.9, -165.3,};
292 //Chamber(5).SetChamberTransform(pos6[0],pos6[1],pos6[2],
293 new TRotMatrix("rot998","rot998",90,0,109.3099,90,19.30999,90);
294 node = new TNode("RICH6","RICH6","S_RICH",pos6[0],pos6[1],pos6[2],"rot998");
297 node->SetLineColor(kColorRICH);
300 Float_t pos7[3]={-161.399,443.3999,-165.3};
301 //Chamber(6).SetChamberTransform(pos7[0],pos7[1],pos7[2],
302 new TRotMatrix("rot999","rot999",90,20,108.1999,110,18.2,110);
303 node = new TNode("RICH7","RICH7","S_RICH",pos7[0],pos7[1],pos7[2],"rot999");
304 node->SetLineColor(kColorRICH);
309 //___________________________________________
310 Int_t AliRICH::DistancetoPrimitive(Int_t , Int_t )
318 //___________________________________________
319 void AliRICH::MakeBranch(Option_t* option)
321 // Create Tree branches for the RICH.
323 const Int_t kBufferSize = 4000;
327 AliDetector::MakeBranch(option);
328 sprintf(branchname,"%sCerenkov",GetName());
329 if (fCerenkovs && gAlice->TreeH()) {
330 gAlice->TreeH()->Branch(branchname,&fCerenkovs, kBufferSize);
331 printf("Making Branch %s for Cerenkov Hits\n",branchname);
334 sprintf(branchname,"%sPadHits",GetName());
335 if (fPadHits && gAlice->TreeH()) {
336 gAlice->TreeH()->Branch(branchname,&fPadHits, kBufferSize);
337 printf("Making Branch %s for PadHits\n",branchname);
340 // one branch for digits per chamber
343 for (i=0; i<kNCH ;i++) {
344 sprintf(branchname,"%sDigits%d",GetName(),i+1);
346 if (fDchambers && gAlice->TreeD()) {
347 gAlice->TreeD()->Branch(branchname,&((*fDchambers)[i]), kBufferSize);
348 printf("Making Branch %s for digits in chamber %d\n",branchname,i+1);
352 // one branch for raw clusters per chamber
353 for (i=0; i<kNCH ;i++) {
354 sprintf(branchname,"%sRawClusters%d",GetName(),i+1);
356 if (fRawClusters && gAlice->TreeR()) {
357 gAlice->TreeR()->Branch(branchname,&((*fRawClusters)[i]), kBufferSize);
358 printf("Making Branch %s for raw clusters in chamber %d\n",branchname,i+1);
362 // one branch for rec hits per chamber
363 for (i=0; i<kNCH ;i++) {
364 sprintf(branchname,"%sRecHits%d",GetName(),i+1);
366 if (fRecHits && gAlice->TreeR()) {
367 gAlice->TreeR()->Branch(branchname,&((*fRecHits)[i]), kBufferSize);
368 printf("Making Branch %s for rec. hits in chamber %d\n",branchname,i+1);
373 //___________________________________________
374 void AliRICH::SetTreeAddress()
376 // Set branch address for the Hits and Digits Tree.
380 AliDetector::SetTreeAddress();
383 TTree *treeH = gAlice->TreeH();
384 TTree *treeD = gAlice->TreeD();
385 TTree *treeR = gAlice->TreeR();
389 branch = treeH->GetBranch("RICHPadHits");
390 if (branch) branch->SetAddress(&fPadHits);
393 branch = treeH->GetBranch("RICHCerenkov");
394 if (branch) branch->SetAddress(&fCerenkovs);
399 for (int i=0; i<kNCH; i++) {
400 sprintf(branchname,"%sDigits%d",GetName(),i+1);
402 branch = treeD->GetBranch(branchname);
403 if (branch) branch->SetAddress(&((*fDchambers)[i]));
408 for (i=0; i<kNCH; i++) {
409 sprintf(branchname,"%sRawClusters%d",GetName(),i+1);
411 branch = treeR->GetBranch(branchname);
412 if (branch) branch->SetAddress(&((*fRawClusters)[i]));
416 for (i=0; i<kNCH; i++) {
417 sprintf(branchname,"%sRecHits%d",GetName(),i+1);
419 branch = treeR->GetBranch(branchname);
420 if (branch) branch->SetAddress(&((*fRecHits)[i]));
426 //___________________________________________
427 void AliRICH::ResetHits()
429 // Reset number of clusters and the cluster array for this detector
430 AliDetector::ResetHits();
433 if (fPadHits) fPadHits->Clear();
434 if (fCerenkovs) fCerenkovs->Clear();
438 //____________________________________________
439 void AliRICH::ResetDigits()
442 // Reset number of digits and the digits array for this detector
444 for ( int i=0;i<kNCH;i++ ) {
445 if ((*fDchambers)[i]) (*fDchambers)[i]->Clear();
446 if (fNdch) fNdch[i]=0;
450 //____________________________________________
451 void AliRICH::ResetRawClusters()
454 // Reset number of raw clusters and the raw clust array for this detector
456 for ( int i=0;i<kNCH;i++ ) {
457 if ((*fRawClusters)[i]) ((TClonesArray*)(*fRawClusters)[i])->Clear();
458 if (fNrawch) fNrawch[i]=0;
462 //____________________________________________
463 void AliRICH::ResetRecHits()
466 // Reset number of raw clusters and the raw clust array for this detector
469 for ( int i=0;i<kNCH;i++ ) {
470 if ((*fRecHits)[i]) ((TClonesArray*)(*fRecHits)[i])->Clear();
471 if (fNrechits) fNrechits[i]=0;
475 //___________________________________________
476 void AliRICH::SetGeometryModel(Int_t id, AliRICHGeometry *geometry)
480 // Setter for the RICH geometry model
484 ((AliRICHChamber*) (*fChambers)[id])->GeometryModel(geometry);
487 //___________________________________________
488 void AliRICH::SetSegmentationModel(Int_t id, AliRICHSegmentation *segmentation)
492 // Setter for the RICH segmentation model
495 ((AliRICHChamber*) (*fChambers)[id])->SegmentationModel(segmentation);
498 //___________________________________________
499 void AliRICH::SetResponseModel(Int_t id, AliRICHResponse *response)
503 // Setter for the RICH response model
506 ((AliRICHChamber*) (*fChambers)[id])->ResponseModel(response);
509 void AliRICH::SetReconstructionModel(Int_t id, AliRICHClusterFinder *reconst)
513 // Setter for the RICH reconstruction model (clusters)
516 ((AliRICHChamber*) (*fChambers)[id])->ReconstructionModel(reconst);
519 void AliRICH::SetNsec(Int_t id, Int_t nsec)
523 // Sets the number of padplanes
526 ((AliRICHChamber*) (*fChambers)[id])->SetNsec(nsec);
530 //___________________________________________
532 void AliRICH::StepManager()
535 // Dummy step manager (should never be called)
539 void AliRICH::FindClusters(Int_t nev,Int_t lastEntry)
543 // Loop on chambers and on cathode planes
545 for (Int_t icat=1;icat<2;icat++) {
546 gAlice->ResetDigits();
547 gAlice->TreeD()->GetEvent(1); // spurious +1 ...
548 for (Int_t ich=0;ich<kNCH;ich++) {
549 AliRICHChamber* iChamber=(AliRICHChamber*) (*fChambers)[ich];
550 TClonesArray *pRICHdigits = this->DigitsAddress(ich);
551 if (pRICHdigits == 0)
554 // Get ready the current chamber stuff
556 AliRICHResponse* response = iChamber->GetResponseModel();
557 AliRICHSegmentation* seg = iChamber->GetSegmentationModel();
558 AliRICHClusterFinder* rec = iChamber->GetReconstructionModel();
560 rec->SetSegmentation(seg);
561 rec->SetResponse(response);
562 rec->SetDigits(pRICHdigits);
563 rec->SetChamber(ich);
564 if (nev==0) rec->CalibrateCOG();
565 rec->FindRawClusters();
568 fRch=RawClustAddress(ich);
572 gAlice->TreeR()->Fill();
574 for (int i=0;i<kNCH;i++) {
575 fRch=RawClustAddress(i);
576 int nraw=fRch->GetEntriesFast();
577 printf ("Chamber %d, raw clusters %d\n",i,nraw);
585 sprintf(hname,"TreeR%d",nev);
586 gAlice->TreeR()->Write(hname,kOverwrite,0);
587 gAlice->TreeR()->Reset();
589 //gObjectTable->Print();
593 //______________________________________________________________________________
594 void AliRICH::Streamer(TBuffer &R__b)
596 // Stream an object of class AliRICH.
597 AliRICHChamber *iChamber;
598 AliRICHSegmentation *segmentation;
599 AliRICHResponse *response;
600 TClonesArray *digitsaddress;
601 TClonesArray *rawcladdress;
602 TClonesArray *rechitaddress;
604 if (R__b.IsReading()) {
605 Version_t R__v = R__b.ReadVersion(); if (R__v) { }
606 AliDetector::Streamer(R__b);
608 R__b >> fPadHits; // diff
610 R__b >> fCerenkovs; // diff
612 R__b >> fRawClusters;
613 R__b >> fRecHits; //diff
614 R__b >> fDebugLevel; //diff
615 R__b.ReadArray(fNdch);
616 R__b.ReadArray(fNrawch);
617 R__b.ReadArray(fNrechits);
620 // Stream chamber related information
621 for (Int_t i =0; i<kNCH; i++) {
622 iChamber=(AliRICHChamber*) (*fChambers)[i];
623 iChamber->Streamer(R__b);
624 segmentation=iChamber->GetSegmentationModel();
625 segmentation->Streamer(R__b);
626 response=iChamber->GetResponseModel();
627 response->Streamer(R__b);
628 rawcladdress=(TClonesArray*) (*fRawClusters)[i];
629 rawcladdress->Streamer(R__b);
630 rechitaddress=(TClonesArray*) (*fRecHits)[i];
631 rechitaddress->Streamer(R__b);
632 digitsaddress=(TClonesArray*) (*fDchambers)[i];
633 digitsaddress->Streamer(R__b);
637 R__b.WriteVersion(AliRICH::IsA());
638 AliDetector::Streamer(R__b);
640 R__b << fPadHits; // diff
642 R__b << fCerenkovs; // diff
644 R__b << fRawClusters;
645 R__b << fRecHits; //diff
646 R__b << fDebugLevel; //diff
647 R__b.WriteArray(fNdch, kNCH);
648 R__b.WriteArray(fNrawch, kNCH);
649 R__b.WriteArray(fNrechits, kNCH);
652 // Stream chamber related information
653 for (Int_t i =0; i<kNCH; i++) {
654 iChamber=(AliRICHChamber*) (*fChambers)[i];
655 iChamber->Streamer(R__b);
656 segmentation=iChamber->GetSegmentationModel();
657 segmentation->Streamer(R__b);
658 response=iChamber->GetResponseModel();
659 response->Streamer(R__b);
660 rawcladdress=(TClonesArray*) (*fRawClusters)[i];
661 rawcladdress->Streamer(R__b);
662 rechitaddress=(TClonesArray*) (*fRecHits)[i];
663 rechitaddress->Streamer(R__b);
664 digitsaddress=(TClonesArray*) (*fDchambers)[i];
665 digitsaddress->Streamer(R__b);
669 AliRICHPadHit* AliRICH::FirstPad(AliRICHHit* hit,TClonesArray *clusters )
672 // Initialise the pad iterator
673 // Return the address of the first padhit for hit
674 TClonesArray *theClusters = clusters;
675 Int_t nclust = theClusters->GetEntriesFast();
676 if (nclust && hit->fPHlast > 0) {
677 sMaxIterPad=Int_t(hit->fPHlast);
678 sCurIterPad=Int_t(hit->fPHfirst);
679 return (AliRICHPadHit*) clusters->UncheckedAt(sCurIterPad-1);
686 AliRICHPadHit* AliRICH::NextPad(TClonesArray *clusters)
689 // Iterates over pads
692 if (sCurIterPad <= sMaxIterPad) {
693 return (AliRICHPadHit*) clusters->UncheckedAt(sCurIterPad-1);
700 void AliRICH::Digitise(Int_t nev, Int_t flag, Option_t *option,Text_t *filename)
702 // keep galice.root for signal and name differently the file for
703 // background when add! otherwise the track info for signal will be lost !
705 static Bool_t first=kTRUE;
707 char *addBackground = strstr(option,"Add");
709 FILE* points; //these will be the digits...
711 points=fopen("points.dat","w");
713 AliRICHChamber* iChamber;
714 AliRICHSegmentation* segmentation;
719 TObjArray *list=new TObjArray;
720 static TClonesArray *pAddress=0;
721 if(!pAddress) pAddress=new TClonesArray("TVector",1000);
724 AliRICH *pRICH = (AliRICH *) gAlice->GetDetector("RICH");
725 AliRICHHitMap* pHitMap[10];
727 for (i=0; i<10; i++) {pHitMap[i]=0;}
728 if (addBackground ) {
731 cout<<"filename"<<fFileName<<endl;
732 pFile=new TFile(fFileName);
733 cout<<"I have opened "<<fFileName<<" file "<<endl;
734 fHits2 = new TClonesArray("AliRICHHit",1000 );
735 fClusters2 = new TClonesArray("AliRICHPadHit",10000);
739 // Get Hits Tree header from file
740 if(fHits2) fHits2->Clear();
741 if(fClusters2) fClusters2->Clear();
742 if(TrH1) delete TrH1;
746 sprintf(treeName,"TreeH%d",nev);
747 TrH1 = (TTree*)gDirectory->Get(treeName);
749 printf("ERROR: cannot find Hits Tree for event:%d\n",nev);
751 // Set branch addresses
754 sprintf(branchname,"%s",GetName());
755 if (TrH1 && fHits2) {
756 branch = TrH1->GetBranch(branchname);
757 if (branch) branch->SetAddress(&fHits2);
759 if (TrH1 && fClusters2) {
760 branch = TrH1->GetBranch("RICHCluster");
761 if (branch) branch->SetAddress(&fClusters2);
768 for (i =0; i<kNCH; i++) {
769 iChamber=(AliRICHChamber*) (*fChambers)[i];
770 segmentation=iChamber->GetSegmentationModel(1);
771 pHitMap[i] = new AliRICHHitMapA1(segmentation, list);
777 TTree *treeH = gAlice->TreeH();
778 Int_t ntracks =(Int_t) treeH->GetEntries();
779 for (Int_t track=0; track<ntracks; track++) {
781 treeH->GetEvent(track);
784 for(AliRICHHit* mHit=(AliRICHHit*)pRICH->FirstHit(-1);
786 mHit=(AliRICHHit*)pRICH->NextHit())
791 Int_t nch = mHit->fChamber-1; // chamber number
792 if (nch >kNCH) continue;
793 iChamber = &(pRICH->Chamber(nch));
795 TParticle *current = (TParticle*)(*gAlice->Particles())[track];
797 Int_t particle = current->GetPdgCode();
799 //printf("Flag:%d\n",flag);
800 //printf("Track:%d\n",track);
801 //printf("Particle:%d\n",particle);
807 if(TMath::Abs(particle) == 211 || TMath::Abs(particle) == 111)
811 if(TMath::Abs(particle)==321 || TMath::Abs(particle)==130 || TMath::Abs(particle)==310
812 || TMath::Abs(particle)==311)
815 if (flag == 3 && TMath::Abs(particle)==2212)
818 if (flag == 4 && TMath::Abs(particle)==13)
821 if (flag == 5 && TMath::Abs(particle)==11)
824 if (flag == 6 && TMath::Abs(particle)==2112)
828 //printf ("Particle: %d, Flag: %d, Digitse: %d\n",particle,flag,digitse);
835 // Loop over pad hits
836 for (AliRICHPadHit* mPad=
837 (AliRICHPadHit*)pRICH->FirstPad(mHit,fPadHits);
839 mPad=(AliRICHPadHit*)pRICH->NextPad(fPadHits))
841 Int_t cathode = mPad->fCathode; // cathode number
842 Int_t ipx = mPad->fPadX; // pad number on X
843 Int_t ipy = mPad->fPadY; // pad number on Y
844 Int_t iqpad = mPad->fQpad; // charge per pad
847 //printf("X:%d, Y:%d, Q:%d\n",ipx,ipy,iqpad);
850 segmentation=iChamber->GetSegmentationModel(cathode);
851 segmentation->GetPadCxy(ipx,ipy,thex,they);
852 new((*pAddress)[countadr++]) TVector(2);
853 TVector &trinfo=*((TVector*) (*pAddress)[countadr-1]);
854 trinfo(0)=(Float_t)track;
855 trinfo(1)=(Float_t)iqpad;
861 AliRICHTransientDigit* pdigit;
862 // build the list of fired pads and update the info
863 if (!pHitMap[nch]->TestHit(ipx, ipy)) {
864 list->AddAtAndExpand(new AliRICHTransientDigit(nch,digits),counter);
865 pHitMap[nch]->SetHit(ipx, ipy, counter);
867 pdigit=(AliRICHTransientDigit*)list->At(list->GetLast());
869 TObjArray *trlist=(TObjArray*)pdigit->TrackList();
870 trlist->Add(&trinfo);
872 pdigit=(AliRICHTransientDigit*) pHitMap[nch]->GetHit(ipx, ipy);
874 (*pdigit).fSignal+=iqpad;
875 // update list of tracks
876 TObjArray* trlist=(TObjArray*)pdigit->TrackList();
877 Int_t lastEntry=trlist->GetLast();
878 TVector *ptrkP=(TVector*)trlist->At(lastEntry);
879 TVector &ptrk=*ptrkP;
880 Int_t lastTrack=Int_t(ptrk(0));
881 Int_t lastCharge=Int_t(ptrk(1));
882 if (lastTrack==track) {
884 trlist->RemoveAt(lastEntry);
886 trinfo(1)=lastCharge;
887 trlist->AddAt(&trinfo,lastEntry);
889 trlist->Add(&trinfo);
891 // check the track list
892 Int_t nptracks=trlist->GetEntriesFast();
894 printf("Attention - tracks: %d (>2)\n",nptracks);
895 //printf("cat,nch,ix,iy %d %d %d %d \n",icat+1,nch,ipx,ipy);
896 for (Int_t tr=0;tr<nptracks;tr++) {
897 TVector *pptrkP=(TVector*)trlist->At(tr);
898 TVector &pptrk=*pptrkP;
899 trk[tr]=Int_t(pptrk(0));
900 chtrk[tr]=Int_t(pptrk(1));
904 } //end loop over clusters
905 }// track type condition
909 // open the file with background
911 if (addBackground ) {
912 ntracks =(Int_t)TrH1->GetEntries();
913 //printf("background - icat,ntracks1 %d %d\n",icat,ntracks);
914 //printf("background - Start loop over tracks \n");
918 for (Int_t trak=0; trak<ntracks; trak++) {
919 if (fHits2) fHits2->Clear();
920 if (fClusters2) fClusters2->Clear();
921 TrH1->GetEvent(trak);
925 for(int j=0;j<fHits2->GetEntriesFast();++j)
927 mHit=(AliRICHHit*) (*fHits2)[j];
928 Int_t nch = mHit->fChamber-1; // chamber number
929 if (nch >6) continue;
930 iChamber = &(pRICH->Chamber(nch));
931 Int_t rmin = (Int_t)iChamber->RInner();
932 Int_t rmax = (Int_t)iChamber->ROuter();
934 // Loop over pad hits
935 for (AliRICHPadHit* mPad=
936 (AliRICHPadHit*)pRICH->FirstPad(mHit,fClusters2);
938 mPad=(AliRICHPadHit*)pRICH->NextPad(fClusters2))
940 Int_t cathode = mPad->fCathode; // cathode number
941 Int_t ipx = mPad->fPadX; // pad number on X
942 Int_t ipy = mPad->fPadY; // pad number on Y
943 Int_t iqpad = mPad->fQpad; // charge per pad
946 segmentation=iChamber->GetSegmentationModel(cathode);
947 segmentation->GetPadCxy(ipx,ipy,thex,they);
948 Float_t rpad=TMath::Sqrt(thex*thex+they*they);
949 if (rpad < rmin || iqpad ==0 || rpad > rmax) continue;
950 new((*pAddress)[countadr++]) TVector(2);
951 TVector &trinfo=*((TVector*) (*pAddress)[countadr-1]);
952 trinfo(0)=-1; // tag background
957 if (trak <4 && nch==0)
958 printf("bgr - pHitMap[nch]->TestHit(ipx, ipy),trak %d %d\n",
959 pHitMap[nch]->TestHit(ipx, ipy),trak);
960 AliRICHTransientDigit* pdigit;
961 // build the list of fired pads and update the info
962 if (!pHitMap[nch]->TestHit(ipx, ipy)) {
963 list->AddAtAndExpand(new AliRICHTransientDigit(nch,digits),counter);
965 pHitMap[nch]->SetHit(ipx, ipy, counter);
967 printf("bgr new elem in list - counter %d\n",counter);
969 pdigit=(AliRICHTransientDigit*)list->At(list->GetLast());
971 TObjArray *trlist=(TObjArray*)pdigit->TrackList();
972 trlist->Add(&trinfo);
974 pdigit=(AliRICHTransientDigit*) pHitMap[nch]->GetHit(ipx, ipy);
976 (*pdigit).fSignal+=iqpad;
977 // update list of tracks
978 TObjArray* trlist=(TObjArray*)pdigit->TrackList();
979 Int_t lastEntry=trlist->GetLast();
980 TVector *ptrkP=(TVector*)trlist->At(lastEntry);
981 TVector &ptrk=*ptrkP;
982 Int_t lastTrack=Int_t(ptrk(0));
986 trlist->Add(&trinfo);
988 // check the track list
989 Int_t nptracks=trlist->GetEntriesFast();
991 for (Int_t tr=0;tr<nptracks;tr++) {
992 TVector *pptrkP=(TVector*)trlist->At(tr);
993 TVector &pptrk=*pptrkP;
994 trk[tr]=Int_t(pptrk(0));
995 chtrk[tr]=Int_t(pptrk(1));
999 } //end loop over clusters
1002 TTree *fAli=gAlice->TreeK();
1003 if (fAli) pFile =fAli->GetCurrentFile();
1009 //cout<<"Start filling digits \n "<<endl;
1010 Int_t nentries=list->GetEntriesFast();
1011 //printf(" \n \n nentries %d \n",nentries);
1013 // start filling the digits
1015 for (Int_t nent=0;nent<nentries;nent++) {
1016 AliRICHTransientDigit *address=(AliRICHTransientDigit*)list->At(nent);
1017 if (address==0) continue;
1019 Int_t ich=address->fChamber;
1020 Int_t q=address->fSignal;
1021 iChamber=(AliRICHChamber*) (*fChambers)[ich];
1022 AliRICHResponse * response=iChamber->GetResponseModel();
1023 Int_t adcmax= (Int_t) response->MaxAdc();
1026 // add white noise and do zero-suppression and signal truncation (new electronics,old electronics gaus 1.2,0.2)
1027 //printf("Treshold: %d\n",iChamber->fTresh->GetHitIndex(address->fPadX,address->fPadY));
1028 Int_t pedestal = iChamber->fTresh->GetHitIndex(address->fPadX,address->fPadY);
1030 //printf("Pedestal:%d\n",pedestal);
1032 Float_t treshold = (pedestal + 4*1.7);
1034 Float_t meanNoise = gRandom->Gaus(1.7, 0.25);
1035 Float_t noise = gRandom->Gaus(0, meanNoise);
1036 q+=(Int_t)(noise + pedestal);
1037 //q+=(Int_t)(noise);
1038 // magic number to be parametrised !!!
1045 if ( q >= adcmax) q=adcmax;
1046 digits[0]=address->fPadX;
1047 digits[1]=address->fPadY;
1050 TObjArray* trlist=(TObjArray*)address->TrackList();
1051 Int_t nptracks=trlist->GetEntriesFast();
1053 // this was changed to accomodate the real number of tracks
1054 if (nptracks > 10) {
1055 cout<<"Attention - tracks > 10 "<<nptracks<<endl;
1059 printf("Attention - tracks > 2 %d \n",nptracks);
1060 //printf("cat,ich,ix,iy,q %d %d %d %d %d \n",
1061 //icat,ich,digits[0],digits[1],q);
1063 for (Int_t tr=0;tr<nptracks;tr++) {
1064 TVector *ppP=(TVector*)trlist->At(tr);
1066 tracks[tr]=Int_t(pp(0));
1067 charges[tr]=Int_t(pp(1));
1068 } //end loop over list of tracks for one pad
1069 if (nptracks < 10 ) {
1070 for (Int_t t=nptracks; t<10; t++) {
1077 fprintf(points,"%4d, %4d, %4d\n",digits[0],digits[1],digits[2]);
1080 pRICH->AddDigits(ich,tracks,charges,digits);
1082 gAlice->TreeD()->Fill();
1085 for(Int_t ii=0;ii<kNCH;++ii) {
1093 //TTree *TD=gAlice->TreeD();
1094 //Stat_t ndig=TD->GetEntries();
1095 //cout<<"number of digits "<<ndig<<endl;
1097 for (int k=0;k<kNCH;k++) {
1098 fDch= pRICH->DigitsAddress(k);
1099 int ndigit=fDch->GetEntriesFast();
1100 printf ("Chamber %d digits %d \n",k,ndigit);
1102 pRICH->ResetDigits();
1104 sprintf(hname,"TreeD%d",nev);
1105 gAlice->TreeD()->Write(hname,kOverwrite,0);
1108 // gAlice->TreeD()->Reset();
1111 // gObjectTable->Print();
1114 AliRICH& AliRICH::operator=(const AliRICH& rhs)
1116 // Assignment operator
1121 Int_t AliRICH::MakePadHits(Float_t xhit,Float_t yhit,Float_t eloss, Int_t idvol, ResponseType res)
1124 // Calls the charge disintegration method of the current chamber and adds
1125 // the simulated cluster to the root treee
1128 Float_t newclust[6][500];
1132 // Integrated pulse height on chamber
1136 ((AliRICHChamber*) (*fChambers)[idvol])->DisIntegration(eloss, xhit, yhit, nnew, newclust, res);
1141 for (Int_t i=0; i<nnew; i++) {
1142 if (Int_t(newclust[3][i]) > 0) {
1145 clhits[1] = Int_t(newclust[5][i]);
1147 clhits[2] = Int_t(newclust[0][i]);
1149 clhits[3] = Int_t(newclust[1][i]);
1151 clhits[4] = Int_t(newclust[2][i]);
1153 clhits[5] = Int_t(newclust[3][i]);
1154 // Pad: chamber sector
1155 clhits[6] = Int_t(newclust[4][i]);