X-Git-Url: http://git.uio.no/git/?a=blobdiff_plain;f=MUON%2FAliMUON.cxx;h=d07d1f94e66ece39772e8cb2929486fd6a7678cc;hb=ecfa008b01fe5050150264fa0daea14d8049e98e;hp=c8c511995e58d81421f0c2fa40fd008ba045d668;hpb=dbe37a89b117331170cb40d2ff7a122c9d431ed4;p=u%2Fmrichter%2FAliRoot.git diff --git a/MUON/AliMUON.cxx b/MUON/AliMUON.cxx index c8c511995e5..d07d1f94e66 100644 --- a/MUON/AliMUON.cxx +++ b/MUON/AliMUON.cxx @@ -12,41 +12,123 @@ * about the suitability of this software for any purpose. It is * * provided "as is" without express or implied warranty. * **************************************************************************/ - /* $Log$ -Revision 1.16 2000/04/26 10:17:31 fca -Changes in Lego for G4 compatibility - -Revision 1.15 2000/01/19 17:16:56 fca -Introducing a list of lists of hits -- more hits allowed for detector now - -Revision 1.14 1999/11/03 13:17:07 fca -Have ProdProcess return const char* - -Revision 1.13 1999/10/26 06:04:48 fca -Introduce TLorentzVector in AliMC::GetSecondary. Thanks to I.Hrivnacova +Revision 1.29 2000/07/28 13:49:38 morsch +SetAcceptance defines inner and outer chamber radii according to angular acceptance. +Can be used for simple acceptance studies. + +Revision 1.28 2000/07/22 16:43:15 morsch +Same comment as before, but now done correctly I hope (sorry it's Saturday evening) + +Revision 1.27 2000/07/22 16:36:50 morsch +Change order of indices in creation (new) of xhit and yhit + +Revision 1.26 2000/07/03 11:54:57 morsch +AliMUONSegmentation and AliMUONHitMap have been replaced by AliSegmentation and AliHitMap in STEER +The methods GetPadIxy and GetPadXxy of AliMUONSegmentation have changed name to GetPadI and GetPadC. + +Revision 1.25 2000/06/29 12:34:09 morsch +AliMUONSegmentation class has been made independent of AliMUONChamber. This makes +it usable with any other geometry class. The link to the object to which it belongs is +established via an index. This assumes that there exists a global geometry manager +from which the pointer to the parent object can be obtained (in our case gAlice). + +Revision 1.24 2000/06/28 15:16:35 morsch +(1) Client code adapted to new method signatures in AliMUONSegmentation (see comments there) +to allow development of slat-muon chamber simulation and reconstruction code in the MUON +framework. The changes should have no side effects (mostly dummy arguments). +(2) Hit disintegration uses 3-dim hit coordinates to allow simulation +of chambers with overlapping modules (MakePadHits, Disintegration). + +Revision 1.23 2000/06/28 12:19:17 morsch +More consequent seperation of global input data services (AliMUONClusterInput singleton) and the +cluster and hit reconstruction algorithms in AliMUONClusterFinderVS. +AliMUONClusterFinderVS becomes the base class for clustering and hit reconstruction. +It requires two cathode planes. Small modifications in the code will make it usable for +one cathode plane and, hence, more general (for test beam data). +AliMUONClusterFinder is now obsolete. + +Revision 1.22 2000/06/28 08:06:10 morsch +Avoid global variables in AliMUONClusterFinderVS by seperating the input data for the fit from the +algorithmic part of the class. Input data resides inside the AliMUONClusterInput singleton. +It also naturally takes care of the TMinuit instance. + +Revision 1.21 2000/06/27 08:54:41 morsch +Problems with on constant array sizes (in hitMap, nmuon, xhit, yhit) corrected. + +Revision 1.20 2000/06/26 14:02:38 morsch +Add class AliMUONConstants with MUON specific constants using static memeber data and access methods. + +Revision 1.19 2000/06/22 13:40:51 morsch +scope problem on HP, "i" declared once +pow changed to TMath::Power (PH, AM) + +Revision 1.18 2000/06/15 07:58:48 morsch +Code from MUON-dev joined + +Revision 1.14.4.17 2000/06/14 14:36:46 morsch +- add TriggerCircuit (PC) +- add GlobalTrigger and LocalTrigger and specific methods (PC) + +Revision 1.14.4.16 2000/06/09 21:20:28 morsch +Most coding rule violations corrected + +Revision 1.14.4.15 2000/05/02 09:54:32 morsch +RULE RN17 violations corrected + +Revision 1.14.4.12 2000/04/26 12:25:02 morsch +Code revised by P. Crochet: +- Z position of TriggerChamber changed according to A.Tournaire Priv.Comm. +- ToF included in the method MakePadHits +- inner radius of flange between beam shielding and trigger corrected +- Trigger global volume updated (according to the new geometry) + +Revision 1.14.4.11 2000/04/19 19:42:08 morsch +Some changes of variable names curing viols and methods concerning +correlated clusters removed. + +Revision 1.14.4.10 2000/03/22 16:44:07 gosset +Memory leak suppressed in function Digitise: +p_adr->Delete() instead of Clear (I.Chevrot and A.Baldisseri) + +Revision 1.14.4.9 2000/03/20 18:15:25 morsch +Positions of trigger chambers corrected (P.C.) + +Revision 1.14.4.8 2000/02/21 15:38:01 morsch +Call to AddHitList introduced to make this version compatible with head. + +Revision 1.14.4.7 2000/02/20 07:45:53 morsch +Bugs in Trigger part of BuildGeomemetry corrected (P.C) + +Revision 1.14.4.6 2000/02/17 14:28:54 morsch +Trigger included into initialization and digitization + +Revision 1.14.4.5 2000/02/15 10:02:58 morsch +Log messages of previous revisions added + +Revision 1.14.4.2 2000/02/04 10:57:34 gosset +Z position of the chambers: +it was the Z position of the stations; +it is now really the Z position of the chambers. + !!!! WARNING: THE CALLS TO "AliMUONChamber::SetZPOS" + !!!! AND "AliMUONChamber::ZPosition" + !!!! HAVE TO BE CHANGED TO "AliMUONChamber::"SetZ" + !!!! AND "AliMUONChamber::Z" -Revision 1.12 1999/10/07 21:08:10 fca -Corrections by G.Chabratova - -Revision 1.11 1999/10/05 17:15:45 fca -Minor syntax for the Alpha OSF - -Revision 1.10 1999/10/01 09:24:40 fca -Protect against no current file in FinishEvent - -Revision 1.9 1999/09/29 09:24:20 fca -Introduction of the Copyright and cvs Log +Revision 1.14.4.3 2000/02/04 16:19:04 gosset +Correction for mis-spelling of NCH + +Revision 1.14.4.4 2000/02/15 09:43:38 morsch +Log message added */ -//////////////////////////////////////////////// + +/////////////////////////////////////////////// // Manager and hits classes for set:MUON // //////////////////////////////////////////////// -#include -#include #include #include #include @@ -66,127 +148,64 @@ Introduction of the Copyright and cvs Log #include #include #include +#include #include "AliMUON.h" -#include "TTUBE.h" -#include "AliMUONClusterFinder.h" +#include "AliMUONHit.h" +#include "AliMUONPadHit.h" +#include "AliMUONDigit.h" +#include "AliMUONTransientDigit.h" +#include "AliMUONRawCluster.h" +#include "AliMUONLocalTrigger.h" +#include "AliMUONGlobalTrigger.h" +#include "AliMUONTriggerCircuit.h" +#include "AliHitMap.h" +#include "AliMUONHitMapA1.h" +#include "AliMUONChamberTrigger.h" +#include "AliMUONConstants.h" +#include "AliMUONClusterFinderVS.h" +#include "AliMUONTriggerDecision.h" #include "AliRun.h" #include "AliMC.h" +#include "AliMUONClusterInput.h" #include "iostream.h" #include "AliCallf77.h" +#include "AliConst.h" + +// Defaults parameters for Z positions of chambers +// taken from values for "stations" in AliMUON::AliMUON +// const Float_t zch[7]={528, 690., 975., 1249., 1449., 1610, 1710.}; +// and from array "dstation" in AliMUONv1::CreateGeometry +// Float_t dstation[5]={20., 20., 20, 20., 20.}; +// for tracking chambers, +// according to (Z1 = zch - dstation) and (Z2 = zch + dstation) +// for the first and second chambers in the station, respectively, +// and from "DTPLANES" in AliMUONv1::CreateGeometry +// const Float_t DTPLANES = 15.; +// for trigger chambers, +// according to (Z1 = zch) and (Z2 = zch + DTPLANES) +// for the first and second chambers in the station, respectively -#ifndef WIN32 -# define reco_init reco_init_ -# define cutpxz cutpxz_ -# define sigmacut sigmacut_ -# define xpreci xpreci_ -# define ypreci ypreci_ -# define reconstmuon reconstmuon_ -# define trackf_read_geant trackf_read_geant_ -# define trackf_read_spoint trackf_read_spoint_ -# define chfill chfill_ -# define chfill2 chfill2_ -# define chf1 chf1_ -# define chfnt chfnt_ -# define hist_create hist_create_ -# define hist_closed hist_closed_ -# define rndm rndm_ -# define fcn fcn_ -# define trackf_fit trackf_fit_ -# define prec_fit prec_fit_ -# define fcnfit fcnfit_ -# define reco_term reco_term_ -#else -# define reco_init RECO_INIT -# define cutpxz CUTPXZ -# define sigmacut SIGMACUT -# define xpreci XPRECI -# define ypreci YPRECI -# define reconstmuon RECONSTMUON -# define trackf_read_geant TRACKF_READ_GEANT -# define trackf_read_spoint TRACKF_READ_SPOINT -# define chfill CHFILL -# define chfill2 CHFILL2 -# define chf1 CHF1 -# define chfnt CHFNT -# define hist_create HIST_CREATE -# define hist_closed HIST_CLOSED -# define rndm RNDM -# define fcn FCN -# define trackf_fit TRACKF_FIT -# define prec_fit PREC_FIT -# define fcnfit FCNFIT -# define reco_term RECO_TERM -#endif - -extern "C" -{ -void type_of_call reco_init(Double_t &, Double_t &, Double_t &); -void type_of_call reco_term(); -void type_of_call cutpxz(Double_t &); -void type_of_call sigmacut(Double_t &); -void type_of_call xpreci(Double_t &); -void type_of_call ypreci(Double_t &); -void type_of_call reconstmuon(Int_t &, Int_t &, Int_t &, Int_t &, Int_t &); -void type_of_call trackf_read_geant(Int_t *, Double_t *, Double_t *, Double_t *, Int_t *, Int_t *, Double_t *, Double_t *, Double_t *, Double_t *,Int_t &, Double_t *, Double_t *, Double_t *, Int_t &, Int_t &, Double_t *, Double_t *, Double_t *, Double_t *); -void type_of_call trackf_read_spoint(Int_t *, Double_t *, Double_t *, Double_t *, Int_t *, Int_t *, Double_t *, Double_t *, Double_t *, Double_t *,Int_t &, Double_t *, Double_t *, Double_t *, Int_t &, Int_t &, Double_t *, Double_t *, Double_t *, Double_t *); -void type_of_call chfill(Int_t &, Float_t &, Float_t &, Float_t &); -void type_of_call chfill2(Int_t &, Float_t &, Float_t &, Float_t &); -void type_of_call chf1(Int_t &, Float_t &, Float_t &); -void type_of_call chfnt(Int_t &, Int_t &, Int_t *, Int_t *, Float_t *, Float_t *, Float_t *, Float_t *, Float_t *, Float_t *, Float_t *, Float_t *); -void type_of_call hist_create(); -void type_of_call hist_closed(); -void type_of_call fcnf(Int_t &, Double_t *, Double_t &, Double_t *, Int_t); -void type_of_call fcn(Int_t &, Double_t *, Double_t &, Double_t *, Int_t &, Int_t &); -void type_of_call trackf_fit(Int_t &, Double_t *, Double_t *, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &); -void type_of_call prec_fit(Double_t &, Double_t &, Double_t &, Double_t &, Double_t&, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &, Double_t &); -void type_of_call fcnfitf(Int_t &, Double_t *, Double_t &, Double_t *, Int_t); -void type_of_call fcnfit(Int_t &, Double_t *, Double_t &, Double_t *, Int_t &, Int_t &); -Float_t type_of_call rndm() {return gRandom->Rndm();} -} - -void fcnfwrap(Int_t &i1, Double_t *d1, Double_t &d2, - Double_t *d3, Int_t i2) -{ - fcnf(i1,d1,d2,d3,i2); -} - -void fcnfitfwrap(Int_t &i1, Double_t *d1, Double_t &d2, - Double_t *d3, Int_t i2) -{ - fcnfitf(i1,d1,d2,d3,i2); -} - - -// Static variables for the pad-hit iterator routines -static Int_t sMaxIterPad=0; -static Int_t sCurIterPad=0; -static TTree *TrH1; -static TTree *TK1; -static TClonesArray *fHits2; //Listof hits for one track only -static TClonesArray *fClusters2; //List of clusters for one track only -static TClonesArray *fParticles2; //List of particles in the Kine tree ClassImp(AliMUON) //___________________________________________ AliMUON::AliMUON() { - fIshunt = 0; - fHits = 0; - fClusters = 0; - fNclusters = 0; - fDchambers = 0; - fNdch = 0; - fRawClusters= 0; - fNrawch = 0; - fCathCorrel= 0; - fNcorch = 0; - fTreeC = 0; - - // modifs perso - fSPxzCut = 0; - fSSigmaCut = 0; - fSXPrec = 0; - fSYPrec = 0; + fIshunt = 0; + fHits = 0; + fPadHits = 0; + fNPadHits = 0; + fDchambers = 0; + fTriggerCircuits = 0; // cp new design of AliMUONTriggerDecision + fNdch = 0; + fRawClusters = 0; + fNrawch = 0; + fGlobalTrigger = 0; + fNLocalTrigger = 0; + fLocalTrigger = 0; + fNLocalTrigger = 0; + fAccMin = 0.; + fAccMax = 0.; + fAccCut = kFALSE; } //___________________________________________ @@ -198,99 +217,151 @@ AliMUON::AliMUON(const char *name, const char *title) */ //End_Html - - fHits = new TClonesArray("AliMUONhit",1000); + + fHits = new TClonesArray("AliMUONHit",1000); gAlice->AddHitList(fHits); - fClusters = new TClonesArray("AliMUONcluster",10000); - fNclusters = 0; + fPadHits = new TClonesArray("AliMUONPadHit",10000); + fNPadHits = 0; fIshunt = 0; - fNdch = new Int_t[10]; + fNdch = new Int_t[AliMUONConstants::NCh()]; - fDchambers = new TObjArray(10); + fDchambers = new TObjArray(AliMUONConstants::NCh()); Int_t i; - for (i=0; i<10 ;i++) { - (*fDchambers)[i] = new TClonesArray("AliMUONdigit",10000); + for (i=0; iSetGid(0); + // Default values for Z of chambers + chamber->SetZ(AliMUONConstants::DefaultChamberZ(ch)); +// + chamber->InitGeo(AliMUONConstants::DefaultChamberZ(ch)); +// Set chamber inner and outer radius to default + chamber->SetRInner(AliMUONConstants::Dmin(st)/2); + chamber->SetROuter(AliMUONConstants::Dmax(st)/2); +// + } // Chamber stCH (0, 1) in + } // Station st (0...) + fMaxStepGas=0.01; + fMaxStepAlu=0.1; + fMaxDestepGas=-1; + fMaxDestepAlu=-1; +// + fMaxIterPad = 0; + fCurIterPad = 0; +// + fAccMin = 0.; + fAccMax = 0.; + fAccCut = kFALSE; + + // cp new design of AliMUONTriggerDecision + fTriggerCircuits = new TObjArray(AliMUONConstants::NTriggerCircuit()); + for (Int_t circ=0; circGetGeometry()->GetNode("alice"); + top=gAlice->GetGeometry()->GetNode("alice"); // MUON // // z-Positions of Chambers - const Float_t cz[5]={511., 686., 971., 1245., 1445.}; -// -// inner diameter - const Float_t dmi[5]={ 35., 47., 67., 86., 100.}; -// -// outer diameter - const Float_t dma[5]={183., 245., 346., 520., 520.}; + const Float_t kCz[7]={511., 686., 971., 1245., 1445., 1600, 1700.}; +// inner diameter (Xlenght for trigger chamber -> active area) + const Float_t kDmin[7]={ 35., 47., 67., 86., 100., 544., 544.}; +// outer diameter (Ylenght for trigger chamber -> active area) + const Float_t kDmax[7]={183., 245., 346., 520., 520., 612., 612.}; TRotMatrix* rot000 = new TRotMatrix("Rot000"," ", 90, 0, 90, 90, 0, 0); TRotMatrix* rot090 = new TRotMatrix("Rot090"," ", 90, 90, 90,180, 0, 0); TRotMatrix* rot180 = new TRotMatrix("Rot180"," ", 90,180, 90,270, 0, 0); TRotMatrix* rot270 = new TRotMatrix("Rot270"," ", 90,270, 90, 0, 0, 0); - - float rmin, rmax, dx, dy, dz, dr, zpos; - float dzc=4.; - char NameChamber[9], NameSense[9], NameFrame[9], NameNode[8]; - for (Int_t i=0; i<5; i++) { + Float_t rmin, rmax, dx, dy, dz, dr, xpos, ypos, zpos; + Float_t dzc1=4.; // tracking chambers + Float_t dzc2=15.; // trigger chambers + Float_t hole=102.; // x-y hole around beam pipe for trig. chambers + Float_t zscale; // scaling parameter trigger chambers + Float_t halfx, halfy; + char nameChamber[9], nameSense[9], nameFrame[9], nameNode[8]; + char nameSense1[9], nameSense2[9]; + for (Int_t i=0; i<7; i++) { for (Int_t j=0; j<2; j++) { Int_t id=2*i+j+1; - if (j==0) { - zpos=cz[i]-dzc; + if (i<5) { // tracking chambers + if (j==0) { + zpos=kCz[i]-dzc1; + } else { + zpos=kCz[i]+dzc1; + } } else { - zpos=cz[i]+dzc; + if (j==0) { + zpos=kCz[i]; + } else { + zpos=kCz[i]+dzc2; + } } - - - sprintf(NameChamber,"C_MUON%d",id); - sprintf(NameSense,"S_MUON%d",id); - sprintf(NameFrame,"F_MUON%d",id); - rmin = dmi[i]/2.-3; - rmax = dma[i]/2.+3; - new TTUBE(NameChamber,"Mother","void",rmin,rmax,0.25,1.); - rmin = dmi[i]/2.; - rmax = dma[i]/2.; - new TTUBE(NameSense,"Sens. region","void",rmin,rmax,0.25, 1.); - dx=(rmax-rmin)/2; - dy=3.; - dz=0.25; - TBRIK* FMUON = new TBRIK(NameFrame,"Frame","void",dx,dy,dz); - Top->cd(); - sprintf(NameNode,"MUON%d",100+id); - Node = new TNode(NameNode,"ChamberNode",NameChamber,0,0,zpos,""); - Node->SetLineColor(kColorMUON); - fNodes->Add(Node); - Node->cd(); - sprintf(NameNode,"MUON%d",200+id); - Node = new TNode(NameNode,"Sens. Region Node",NameSense,0,0,0,""); - Node->SetLineColor(kColorMUON); - fNodes->Add(Node); - Node->cd(); - dr=dx+rmin; - sprintf(NameNode,"MUON%d",300+id); - NodeF = new TNode(NameNode,"Frame0",FMUON,dr, 0, 0,rot000,""); - NodeF->SetLineColor(kColorMUON); - fNodes->Add(NodeF); - Node->cd(); - sprintf(NameNode,"MUON%d",400+id); - NodeF = new TNode(NameNode,"Frame1",FMUON,0 ,dr,0,rot090,""); - NodeF->SetLineColor(kColorMUON); - fNodes->Add(NodeF); - Node->cd(); - sprintf(NameNode,"MUON%d",500+id); - NodeF = new TNode(NameNode,"Frame2",FMUON,-dr,0,0,rot180,""); - NodeF->SetLineColor(kColorMUON); - fNodes->Add(NodeF); - Node ->cd(); - sprintf(NameNode,"MUON%d",600+id); - NodeF = new TNode(NameNode,"Frame3",FMUON,0,-dr,0,rot270,""); - NodeF->SetLineColor(kColorMUON); - fNodes->Add(NodeF); + sprintf(nameChamber,"C_MUON%d",id); + sprintf(nameSense,"S_MUON%d",id); + sprintf(nameSense1,"S1_MUON%d",id); + sprintf(nameSense2,"S2_MUON%d",id); + sprintf(nameFrame,"F_MUON%d",id); + if (i<5) { // tracking chambers + rmin = kDmin[i]/2.-3; + rmax = kDmax[i]/2.+3; + new TTUBE(nameChamber,"Mother","void",rmin,rmax,0.25,1.); + rmin = kDmin[i]/2.; + rmax = kDmax[i]/2.; + new TTUBE(nameSense,"Sens. region","void",rmin,rmax,0.25, 1.); + dx=(rmax-rmin)/2; + dy=3.; + dz=0.25; + TBRIK* frMUON = new TBRIK(nameFrame,"Frame","void",dx,dy,dz); + top->cd(); + sprintf(nameNode,"MUON%d",100+id); + node = new TNode(nameNode,"ChamberNode",nameChamber,0,0,zpos,""); + node->SetLineColor(kColorMUON); + fNodes->Add(node); + node->cd(); + sprintf(nameNode,"MUON%d",200+id); + node = new TNode(nameNode,"Sens. Region Node",nameSense,0,0,0,""); + node->SetLineColor(kColorMUON); + node->cd(); + dr=dx+rmin; + sprintf(nameNode,"MUON%d",300+id); + nodeF = new TNode(nameNode,"Frame0",frMUON,dr, 0, 0,rot000,""); + nodeF->SetLineColor(kColorMUON); + node->cd(); + sprintf(nameNode,"MUON%d",400+id); + nodeF = new TNode(nameNode,"Frame1",frMUON,0 ,dr,0,rot090,""); + nodeF->SetLineColor(kColorMUON); + node->cd(); + sprintf(nameNode,"MUON%d",500+id); + nodeF = new TNode(nameNode,"Frame2",frMUON,-dr,0,0,rot180,""); + nodeF->SetLineColor(kColorMUON); + node ->cd(); + sprintf(nameNode,"MUON%d",600+id); + nodeF = new TNode(nameNode,"Frame3",frMUON,0,-dr,0,rot270,""); + nodeF->SetLineColor(kColorMUON); + } else { + zscale=zpos/kCz[5]; + Float_t xsize=kDmin[i]*zscale; + Float_t ysize=kDmax[i]*zscale; + Float_t holeScaled=hole*zscale; + + halfx=xsize/2.+3.; + halfy=ysize/2.+3.; + new TBRIK(nameChamber,"Mother","void",halfx,halfy,0.25); + top->cd(); + sprintf(nameNode,"MUON%d",100+id); + node = new TNode(nameNode,"Chambernode",nameChamber,0,0,zpos,""); + node->SetLineColor(kColorMUON2); + fNodes->Add(node); + +// up/down of beam pipe + halfx=xsize/2.; + halfy=(ysize/2.-holeScaled/2.)/2.; + new TBRIK(nameSense,"Sens. region","void",halfx,halfy,0.25); + + node->cd(); + ypos=holeScaled/2.+((ysize/2.-holeScaled/2.)/2.); + sprintf(nameNode,"MUON%d",200+id); + nodeS = new TNode(nameNode,"Sens. Region Node",nameSense,0,ypos,0,""); + nodeS->SetLineColor(kColorMUON2); + + node->cd(); + ypos=-1.*ypos; + sprintf(nameNode,"MUON%d",300+id); + nodeS = new TNode(nameNode,"Sens. Region Node",nameSense,0,ypos,0,""); + nodeS->SetLineColor(kColorMUON2); + +// left/right of beam pipe + halfx=(xsize/2.-holeScaled/2.)/2.; + halfy=holeScaled/2.; + new TBRIK(nameSense1,"Sens. region","void",halfx,halfy,0.25); + + node->cd(); + xpos=holeScaled/2.+((xsize/2.-holeScaled/2.)/2.); + sprintf(nameNode,"MUON%d",400+id); + nodeS = new TNode(nameNode,"Sens. Region Node",nameSense1,xpos,0,0,""); + nodeS->SetLineColor(kColorMUON2); + + node->cd(); + xpos=-1.*xpos; + sprintf(nameNode,"MUON%d",500+id); + nodeS = new TNode(nameNode,"Sens. Region Node",nameSense1,xpos,0,0,""); + nodeS->SetLineColor(kColorMUON2); + +// missing corners + halfx=17.*zscale/2.; + halfy=halfx; + new TBRIK(nameSense2,"Sens. region","void",halfx,halfy,0.25); + + node->cd(); + xpos=holeScaled/2.-halfx; + ypos=xpos; + sprintf(nameNode,"MUON%d",600+id); + nodeS = new TNode(nameNode,"Sens. Region Node",nameSense2,xpos,ypos,0,""); + nodeS->SetLineColor(kColorMUON2); + + node->cd(); + ypos=-1.*xpos; + sprintf(nameNode,"MUON%d",700+id); + nodeS = new TNode(nameNode,"Sens. Region Node",nameSense2,xpos,ypos,0,""); + nodeS->SetLineColor(kColorMUON2); + + node->cd(); + xpos=-1.*xpos; + sprintf(nameNode,"MUON%d",800+id); + nodeS = new TNode(nameNode,"Sens. Region Node",nameSense2,xpos,ypos,0,""); + nodeS->SetLineColor(kColorMUON2); + + node->cd(); + ypos=-1.*xpos; + sprintf(nameNode,"MUON%d",900+id); + nodeS = new TNode(nameNode,"Sens. Region Node",nameSense2,xpos,ypos,0,""); + nodeS->SetLineColor(kColorMUON2); + } } } } @@ -419,43 +585,55 @@ Int_t AliMUON::DistancetoPrimitive(Int_t , Int_t ) //___________________________________________ void AliMUON::MakeBranch(Option_t* option) { - // Create Tree branches for the MUON. - - const Int_t buffersize = 4000; - char branchname[30]; - sprintf(branchname,"%sCluster",GetName()); - - AliDetector::MakeBranch(option); - - if (fClusters && gAlice->TreeH()) { - gAlice->TreeH()->Branch(branchname,&fClusters, buffersize); - printf("Making Branch %s for clusters\n",branchname); - } - + // Create Tree branches for the MUON. + const Int_t kBufferSize = 4000; + char branchname[30]; + sprintf(branchname,"%sCluster",GetName()); + + AliDetector::MakeBranch(option); + + if (fPadHits && gAlice->TreeH()) { + gAlice->TreeH()->Branch(branchname,&fPadHits, kBufferSize); + printf("Making Branch %s for clusters\n",branchname); + } + // one branch for digits per chamber - Int_t i; - - for (i=0; i<10 ;i++) { - sprintf(branchname,"%sDigits%d",GetName(),i+1); - - if (fDchambers && gAlice->TreeD()) { - gAlice->TreeD()->Branch(branchname,&((*fDchambers)[i]), buffersize); - printf("Making Branch %s for digits in chamber %d\n",branchname,i+1); - } - } - - //printf("Make Branch - TreeR address %p\n",gAlice->TreeR()); - + Int_t i; + + for (i=0; iTreeD()) { + gAlice->TreeD()->Branch(branchname,&((*fDchambers)[i]), kBufferSize); + printf("Making Branch %s for digits in chamber %d\n",branchname,i+1); + } + } + + printf("Make Branch - TreeR address %p\n",gAlice->TreeR()); + // one branch for raw clusters per chamber - for (i=0; i<10 ;i++) { - sprintf(branchname,"%sRawClusters%d",GetName(),i+1); - - if (fRawClusters && gAlice->TreeR()) { - gAlice->TreeR()->Branch(branchname,&((*fRawClusters)[i]), buffersize); - printf("Making Branch %s for raw clusters in chamber %d\n",branchname,i+1); - } - } - + for (i=0; iTreeR()) { + gAlice->TreeR()->Branch(branchname,&((*fRawClusters)[i]), kBufferSize); + printf("Making Branch %s for raw clusters in chamber %d\n",branchname,i+1); + } + } + +// one branch for global trigger + sprintf(branchname,"%sGlobalTrigger",GetName()); + if (fGlobalTrigger && gAlice->TreeR()) { + gAlice->TreeR()->Branch(branchname,&fGlobalTrigger,kBufferSize); + printf("Making Branch %s for Global Trigger\n",branchname); + } +// one branch for local trigger + sprintf(branchname,"%sLocalTrigger",GetName()); + if (fLocalTrigger && gAlice->TreeR()) { + gAlice->TreeR()->Branch(branchname,&fLocalTrigger,kBufferSize); + printf("Making Branch %s for Local Trigger\n",branchname); + } + } //___________________________________________ @@ -471,14 +649,14 @@ void AliMUON::SetTreeAddress() TTree *treeR = gAlice->TreeR(); if (treeH) { - if (fClusters) { + if (fPadHits) { branch = treeH->GetBranch("MUONCluster"); - if (branch) branch->SetAddress(&fClusters); + if (branch) branch->SetAddress(&fPadHits); } } if (treeD) { - for (int i=0; i<10; i++) { + for (int i=0; iGetBranch(branchname); @@ -490,23 +668,31 @@ void AliMUON::SetTreeAddress() // printf("SetTreeAddress --- treeR address %p \n",treeR); if (treeR) { - for (int i=0; i<10; i++) { + for (int i=0; iGetBranch(branchname); if (branch) branch->SetAddress(&((*fRawClusters)[i])); } } - } + if (fLocalTrigger) { + branch = treeR->GetBranch("MUONLocalTrigger"); + if (branch) branch->SetAddress(&fLocalTrigger); + } + if (fGlobalTrigger) { + branch = treeR->GetBranch("MUONGlobalTrigger"); + if (branch) branch->SetAddress(&fGlobalTrigger); + } + } } //___________________________________________ void AliMUON::ResetHits() { // Reset number of clusters and the cluster array for this detector AliDetector::ResetHits(); - fNclusters = 0; - if (fClusters) fClusters->Clear(); + fNPadHits = 0; + if (fPadHits) fPadHits->Clear(); } //____________________________________________ @@ -515,7 +701,7 @@ void AliMUON::ResetDigits() // // Reset number of digits and the digits array for this detector // - for ( int i=0;i<10;i++ ) { + for ( int i=0;iClear(); if (fNdch) fNdch[i]=0; } @@ -526,63 +712,79 @@ void AliMUON::ResetRawClusters() // // Reset number of raw clusters and the raw clust array for this detector // - for ( int i=0;i<10;i++ ) { + for ( int i=0;iClear(); if (fNrawch) fNrawch[i]=0; } } + //____________________________________________ -void AliMUON::ResetCorrelation() +void AliMUON::ResetTrigger() { - // - // Reset number of correl clusters and the correl clust array for - // this detector - // - for ( int i=0;i<10;i++ ) { - if ((*fCathCorrel)[i]) ((TClonesArray*)(*fCathCorrel)[i])->Clear(); - if (fNcorch) fNcorch[i]=0; - } + // Reset Local and Global Trigger + fNGlobalTrigger = 0; + if (fGlobalTrigger) fGlobalTrigger->Clear(); + fNLocalTrigger = 0; + if (fLocalTrigger) fLocalTrigger->Clear(); +} + +//____________________________________________ +void AliMUON::SetPadSize(Int_t id, Int_t isec, Float_t p1, Float_t p2) +{ + Int_t i=2*(id-1); + ((AliMUONChamber*) (*fChambers)[i]) ->SetPadSize(isec,p1,p2); + ((AliMUONChamber*) (*fChambers)[i+1])->SetPadSize(isec,p1,p2); } //___________________________________________ +void AliMUON::SetChambersZ(const Float_t *Z) +{ + // Set Z values for all chambers (tracking and trigger) + // from the array pointed to by "Z" + for (Int_t ch = 0; ch < AliMUONConstants::NCh(); ch++) + ((AliMUONChamber*) ((*fChambers)[ch]))->SetZ(Z[ch]); + return; +} -void AliMUON::SetPADSIZ(Int_t id, Int_t isec, Float_t p1, Float_t p2) +//___________________________________________ +void AliMUON::SetChambersZToDefault() { - Int_t i=2*(id-1); - ((AliMUONchamber*) (*fChambers)[i]) ->SetPADSIZ(isec,p1,p2); - ((AliMUONchamber*) (*fChambers)[i+1])->SetPADSIZ(isec,p1,p2); + // Set Z values for all chambers (tracking and trigger) + // to default values + SetChambersZ(AliMUONConstants::DefaultChamberZ()); + return; } //___________________________________________ void AliMUON::SetChargeSlope(Int_t id, Float_t p1) { Int_t i=2*(id-1); - ((AliMUONchamber*) (*fChambers)[i])->SetChargeSlope(p1); - ((AliMUONchamber*) (*fChambers)[i+1])->SetChargeSlope(p1); + ((AliMUONChamber*) (*fChambers)[i])->SetChargeSlope(p1); + ((AliMUONChamber*) (*fChambers)[i+1])->SetChargeSlope(p1); } //___________________________________________ void AliMUON::SetChargeSpread(Int_t id, Float_t p1, Float_t p2) { Int_t i=2*(id-1); - ((AliMUONchamber*) (*fChambers)[i])->SetChargeSpread(p1,p2); - ((AliMUONchamber*) (*fChambers)[i+1])->SetChargeSpread(p1,p2); + ((AliMUONChamber*) (*fChambers)[i])->SetChargeSpread(p1,p2); + ((AliMUONChamber*) (*fChambers)[i+1])->SetChargeSpread(p1,p2); } //___________________________________________ void AliMUON::SetSigmaIntegration(Int_t id, Float_t p1) { Int_t i=2*(id-1); - ((AliMUONchamber*) (*fChambers)[i])->SetSigmaIntegration(p1); - ((AliMUONchamber*) (*fChambers)[i+1])->SetSigmaIntegration(p1); + ((AliMUONChamber*) (*fChambers)[i])->SetSigmaIntegration(p1); + ((AliMUONChamber*) (*fChambers)[i+1])->SetSigmaIntegration(p1); } //___________________________________________ -void AliMUON::SetMaxAdc(Int_t id, Float_t p1) +void AliMUON::SetMaxAdc(Int_t id, Int_t p1) { Int_t i=2*(id-1); - ((AliMUONchamber*) (*fChambers)[i])->SetMaxAdc(p1); - ((AliMUONchamber*) (*fChambers)[i+1])->SetMaxAdc(p1); + ((AliMUONChamber*) (*fChambers)[i])->SetMaxAdc(p1); + ((AliMUONChamber*) (*fChambers)[i+1])->SetMaxAdc(p1); } //___________________________________________ @@ -609,74 +811,53 @@ void AliMUON::SetMaxDestepAlu(Float_t p1) fMaxDestepAlu=p1; } //___________________________________________ -void AliMUON::SetMuonAcc(Bool_t acc, Float_t angmin, Float_t angmax) +void AliMUON::SetAcceptance(Bool_t acc, Float_t angmin, Float_t angmax) { fAccCut=acc; - fAccMin=angmin; - fAccMax=angmax; + fAccMin=angmin*TMath::Pi()/180; + fAccMax=angmax*TMath::Pi()/180; + Int_t ch; + if (acc) { + for (Int_t st = 0; st < AliMUONConstants::NCh() / 2; st++) { + // Loop over 2 chambers in the station + for (Int_t stCH = 0; stCH < 2; stCH++) { + ch = 2 * st + stCH; +// Set chamber inner and outer radius according to acceptance cuts + Chamber(ch).SetRInner(AliMUONConstants::DefaultChamberZ(ch)*TMath::Tan(fAccMin)); + Chamber(ch).SetROuter(AliMUONConstants::DefaultChamberZ(ch)*TMath::Tan(fAccMax)); + } // chamber loop + } // station loop + } } //___________________________________________ -void AliMUON::SetSegmentationModel(Int_t id, Int_t isec, AliMUONsegmentation *segmentation) +void AliMUON::SetSegmentationModel(Int_t id, Int_t isec, AliSegmentation *segmentation) { - ((AliMUONchamber*) (*fChambers)[id])->SegmentationModel(isec, segmentation); + ((AliMUONChamber*) (*fChambers)[id])->SetSegmentationModel(isec, segmentation); } //___________________________________________ -void AliMUON::SetResponseModel(Int_t id, AliMUONresponse *response) +void AliMUON::SetResponseModel(Int_t id, AliMUONResponse *response) { - ((AliMUONchamber*) (*fChambers)[id])->ResponseModel(response); + ((AliMUONChamber*) (*fChambers)[id])->SetResponseModel(response); } -void AliMUON::SetReconstructionModel(Int_t id, AliMUONClusterFinder *reconst) +void AliMUON::SetReconstructionModel(Int_t id, AliMUONClusterFinderVS *reconst) { - ((AliMUONchamber*) (*fChambers)[id])->ReconstructionModel(reconst); + ((AliMUONChamber*) (*fChambers)[id])->SetReconstructionModel(reconst); } void AliMUON::SetNsec(Int_t id, Int_t nsec) { - ((AliMUONchamber*) (*fChambers)[id])->SetNsec(nsec); + ((AliMUONChamber*) (*fChambers)[id])->SetNsec(nsec); } //___________________________________________ -void AliMUON::StepManager() -{ - printf("Dummy version of muon step -- it should never happen!!\n"); - /* - const Float_t kRaddeg = 180/TMath::Pi(); - Int_t nsec, ipart; - TLorentzVector x, p; - Float_t pt, th0, th2; - char *proc; - if(fAccCut) { - if((nsec=gMC->NSecondaries())>0) { - proc=gMC->ProdProcess(); - if((gMC->TrackPid()==443 || gMC->TrackPid()==553) && !strcmp(proc,"DCAY")) { - // - // Check angular acceptance - // --- and have muons from resonance decays in the wanted window --- - if(nsec != 2) { - printf(" AliMUON::StepManager: Strange resonance Decay into %d particles\n",nsec); - gMC->StopEvent(); - } else { - gMC->GetSecondary(0,ipart,x,p); - pt = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]); - th0 = TMath::ATan2(pt,p[2])*kRaddeg; - gMC->GetSecondary(1,ipart,x,p); - pt = TMath::Sqrt(p[0]*p[0]+p[1]*p[1]); - th2 = TMath::ATan2(pt,p[2])*kRaddeg; - if(!(fAccMin < th0 && th0 < fAccMax) || - !(fAccMin < th2 && th2 < fAccMax)) - gMC->StopEvent(); - } - } - } - } - */ -} -void AliMUON::MakePadHits(Float_t xhit,Float_t yhit,Float_t eloss, Int_t idvol) + +void AliMUON::MakePadHits(Float_t xhit,Float_t yhit, Float_t zhit, + Float_t eloss, Float_t tof, Int_t idvol) { // // Calls the charge disintegration method of the current chamber and adds @@ -694,8 +875,8 @@ void AliMUON::MakePadHits(Float_t xhit,Float_t yhit,Float_t eloss, Int_t idvol) clhits[0]=fNhits+1; // // - ((AliMUONchamber*) (*fChambers)[idvol])->DisIntegration(eloss, xhit, yhit, nnew, newclust); -// printf("\n Add new clusters %d %f \n", nnew, eloss*1.e9); + ((AliMUONChamber*) (*fChambers)[idvol]) + ->DisIntegration(eloss, tof, xhit, yhit, zhit, nnew, newclust); Int_t ic=0; // @@ -716,162 +897,162 @@ void AliMUON::MakePadHits(Float_t xhit,Float_t yhit,Float_t eloss, Int_t idvol) // Pad: chamber sector clhits[6] = Int_t(newclust[4][i]); - AddCluster(clhits); + AddPadHit(clhits); } } -// printf("\n %d new clusters added", ic); } -void AliMUON::Digitise(Int_t nev,Int_t bgr_ev,Option_t *option, Option_t *,Text_t *filename) +//---------------------------------------------------------------------- + +void AliMUON::Digitise(Int_t nev,Int_t bgrEvent,Option_t *option,Option_t *opt,Text_t *filename) { // keep galice.root for signal and name differently the file for // background when add! otherwise the track info for signal will be lost ! static Bool_t first=kTRUE; -// static TTree *TrH1; - static TFile *File; - char *Add = strstr(option,"Add"); - //char *listoftracks = strstr(opt,"listoftracks"); + static TFile *file; + char *addBackground = strstr(option,"Add"); - AliMUONchamber* iChamber; - AliMUONsegmentation* segmentation; + + AliMUONChamber* iChamber; + AliSegmentation* segmentation; Int_t trk[50]; Int_t chtrk[50]; TObjArray *list=new TObjArray; - static TClonesArray *p_adr=0; - if(!p_adr) p_adr=new TClonesArray("TVector",1000); + static TClonesArray *pAddress=0; + if(!pAddress) pAddress=new TClonesArray("TVector",1000); Int_t digits[5]; - AliMUON *MUON = (AliMUON *) gAlice->GetModule("MUON"); - AliMUONHitMap * HitMap[10]; - for (Int_t i=0; i<10; i++) {HitMap[i]=0;} - if (Add ) { + AliMUON *pMUON = (AliMUON *) gAlice->GetModule("MUON"); + AliHitMap** hitMap= new AliHitMap* [AliMUONConstants::NCh()]; + for (Int_t i=0; icd(); - //File->ls(); + file->cd(); + //file->ls(); // Get Hits Tree header from file if(fHits2) fHits2->Clear(); - if(fClusters2) fClusters2->Clear(); - if(TrH1) delete TrH1; - TrH1=0; + if(fPadHits2) fPadHits2->Clear(); + if(fTrH1) delete fTrH1; + fTrH1=0; char treeName[20]; - sprintf(treeName,"TreeH%d",bgr_ev); - TrH1 = (TTree*)gDirectory->Get(treeName); - //printf("TrH1 %p of treename %s for event %d \n",TrH1,treeName,bgr_ev); + sprintf(treeName,"TreeH%d",bgrEvent); + fTrH1 = (TTree*)gDirectory->Get(treeName); + //printf("fTrH1 %p of treename %s for event %d \n",fTrH1,treeName,bgrEvent); - if (!TrH1) { - printf("ERROR: cannot find Hits Tree for event:%d\n",bgr_ev); + if (!fTrH1) { + printf("ERROR: cannot find Hits Tree for event:%d\n",bgrEvent); } // Set branch addresses TBranch *branch; char branchname[20]; sprintf(branchname,"%s",GetName()); - if (TrH1 && fHits2) { - branch = TrH1->GetBranch(branchname); + if (fTrH1 && fHits2) { + branch = fTrH1->GetBranch(branchname); if (branch) branch->SetAddress(&fHits2); } - if (TrH1 && fClusters2) { - branch = TrH1->GetBranch("MUONCluster"); - if (branch) branch->SetAddress(&fClusters2); + if (fTrH1 && fPadHits2) { + branch = fTrH1->GetBranch("MUONCluster"); + if (branch) branch->SetAddress(&fPadHits2); } // test - //Int_t ntracks1 =(Int_t)TrH1->GetEntries(); + //Int_t ntracks1 =(Int_t)fTrH1->GetEntries(); //printf("background - ntracks1 - %d\n",ntracks1); } // // loop over cathodes // - AliMUONHitMap* hm; + AliHitMap* hm; Int_t countadr=0; for (int icat=0; icat<2; icat++) { Int_t counter=0; - for (Int_t i =0; i<10; i++) { - iChamber=(AliMUONchamber*) (*fChambers)[i]; + Int_t * nmuon = new Int_t [AliMUONConstants::NCh()]; + for (Int_t i =0; iNsec()==1 && icat==1) { continue; } else { - segmentation=iChamber->GetSegmentationModel(icat+1); + segmentation=iChamber->SegmentationModel(icat+1); } - HitMap[i] = new AliMUONHitMapA1(segmentation, list); + hitMap[i] = new AliMUONHitMapA1(segmentation, list); + nmuon[i]=0; } //printf("Start loop over tracks \n"); // // Loop over tracks // - TTree *TH = gAlice->TreeH(); - Int_t ntracks =(Int_t) TH->GetEntries(); - //printf("signal - ntracks %d\n",ntracks); - Int_t nmuon[10]={0,0,0,0,0,0,0,0,0,0}; - Float_t xhit[10][2]; - Float_t yhit[10][2]; - + TTree *treeH = gAlice->TreeH(); + Int_t ntracks =(Int_t) treeH->GetEntries(); + Int_t jj; + + Float_t ** xhit = new Float_t * [AliMUONConstants::NCh()]; + for (jj=0; jjResetHits(); - TH->GetEvent(track); - + treeH->GetEvent(track); // // Loop over hits - for(AliMUONhit* mHit=(AliMUONhit*)MUON->FirstHit(-1); + for(AliMUONHit* mHit=(AliMUONHit*)pMUON->FirstHit(-1); mHit; - mHit=(AliMUONhit*)MUON->NextHit()) + mHit=(AliMUONHit*)pMUON->NextHit()) { Int_t nch = mHit->fChamber-1; // chamber number - if (nch >9) continue; - iChamber = &(MUON->Chamber(nch)); - Int_t rmin = (Int_t)iChamber->RInner(); - Int_t rmax = (Int_t)iChamber->ROuter(); + if (nch > AliMUONConstants::NCh()-1) continue; + iChamber = &(pMUON->Chamber(nch)); // new 17.07.99 - if (Add) { - - if (mHit->fParticle == kMuonPlus || mHit->fParticle == kMuonMinus) { - xhit[nch][nmuon[nch]]=mHit->fX; - yhit[nch][nmuon[nch]]=mHit->fY; - nmuon[nch]++; - if (nmuon[nch] >2) printf("nmuon %d\n",nmuon[nch]); - - } + if (addBackground) { + + if (mHit->fParticle == kMuonPlus + || mHit->fParticle == kMuonMinus) { + xhit[nch][nmuon[nch]]=mHit->fX; + yhit[nch][nmuon[nch]]=mHit->fY; + nmuon[nch]++; + if (nmuon[nch] >2) printf("nmuon %d\n",nmuon[nch]); + } } - + // // Loop over pad hits - for (AliMUONcluster* mPad= - (AliMUONcluster*)MUON->FirstPad(mHit,fClusters); + for (AliMUONPadHit* mPad= + (AliMUONPadHit*)pMUON->FirstPad(mHit,fPadHits); mPad; - mPad=(AliMUONcluster*)MUON->NextPad(fClusters)) + mPad=(AliMUONPadHit*)pMUON->NextPad(fPadHits)) { Int_t cathode = mPad->fCathode; // cathode number Int_t ipx = mPad->fPadX; // pad number on X Int_t ipy = mPad->fPadY; // pad number on Y - Int_t iqpad = Int_t(mPad->fQpad*kScale);// charge per pad -// Int_t iqpad = mPad->fQpad; // charge per pad + Int_t iqpad = Int_t(mPad->fQpad);// charge per pad // // if (cathode != (icat+1)) continue; // fill the info array - Float_t thex, they; - segmentation=iChamber->GetSegmentationModel(cathode); - segmentation->GetPadCxy(ipx,ipy,thex,they); - Float_t rpad=TMath::Sqrt(thex*thex+they*they); - if (rpad < rmin || iqpad ==0 || rpad > rmax) continue; - - new((*p_adr)[countadr++]) TVector(2); - TVector &trinfo=*((TVector*) (*p_adr)[countadr-1]); + Float_t thex, they, thez; + segmentation=iChamber->SegmentationModel(cathode); + segmentation->GetPadC(ipx,ipy,thex,they,thez); +// Float_t rpad=TMath::Sqrt(thex*thex+they*they); +// if (rpad < rmin || iqpad ==0 || rpad > rmax) continue; + + new((*pAddress)[countadr++]) TVector(2); + TVector &trinfo=*((TVector*) (*pAddress)[countadr-1]); trinfo(0)=(Float_t)track; trinfo(1)=(Float_t)iqpad; @@ -879,41 +1060,42 @@ void AliMUON::Digitise(Int_t nev,Int_t bgr_ev,Option_t *option, Option_t *,Text_ digits[1]=ipy; digits[2]=iqpad; digits[3]=iqpad; - if (mHit->fParticle == kMuonPlus || mHit->fParticle == kMuonMinus) { - digits[4]=mPad->fHitNumber; + if (mHit->fParticle == kMuonPlus || + mHit->fParticle == kMuonMinus) { + digits[4]=mPad->fHitNumber; } else digits[4]=-1; - AliMUONlist* pdigit; + AliMUONTransientDigit* pdigit; // build the list of fired pads and update the info - if (!HitMap[nch]->TestHit(ipx, ipy)) { + if (!hitMap[nch]->TestHit(ipx, ipy)) { list->AddAtAndExpand( - new AliMUONlist(nch,digits),counter); + new AliMUONTransientDigit(nch,digits),counter); - HitMap[nch]->SetHit(ipx, ipy, counter); + hitMap[nch]->SetHit(ipx, ipy, counter); counter++; - pdigit=(AliMUONlist*)list->At(list->GetLast()); + pdigit=(AliMUONTransientDigit*)list->At(list->GetLast()); // list of tracks TObjArray *trlist=(TObjArray*)pdigit->TrackList(); trlist->Add(&trinfo); } else { - pdigit=(AliMUONlist*) HitMap[nch]->GetHit(ipx, ipy); + pdigit=(AliMUONTransientDigit*) hitMap[nch]->GetHit(ipx, ipy); // update charge (*pdigit).fSignal+=iqpad; (*pdigit).fPhysics+=iqpad; // update list of tracks TObjArray* trlist=(TObjArray*)pdigit->TrackList(); - Int_t last_entry=trlist->GetLast(); - TVector *ptrk_p=(TVector*)trlist->At(last_entry); - TVector &ptrk=*ptrk_p; - Int_t last_track=Int_t(ptrk(0)); - Int_t last_charge=Int_t(ptrk(1)); - if (last_track==track) { - last_charge+=iqpad; - trlist->RemoveAt(last_entry); - trinfo(0)=last_track; - trinfo(1)=last_charge; - trlist->AddAt(&trinfo,last_entry); + Int_t lastEntry=trlist->GetLast(); + TVector *pTrack=(TVector*)trlist->At(lastEntry); + TVector &ptrk=*pTrack; + Int_t lastTrack=Int_t(ptrk(0)); + Int_t lastCharge=Int_t(ptrk(1)); + if (lastTrack==track) { + lastCharge+=iqpad; + trlist->RemoveAt(lastEntry); + trinfo(0)=lastTrack; + trinfo(1)=lastCharge; + trlist->AddAt(&trinfo,lastEntry); } else { trlist->Add(&trinfo); } @@ -921,8 +1103,8 @@ void AliMUON::Digitise(Int_t nev,Int_t bgr_ev,Option_t *option, Option_t *,Text_ Int_t nptracks=trlist->GetEntriesFast(); if (nptracks > 2) { for (Int_t tr=0;trAt(tr); - TVector &pptrk=*pptrk_p; + TVector *ppTrack=(TVector*)trlist->At(tr); + TVector &pptrk=*ppTrack; trk[tr]=Int_t(pptrk(0)); chtrk[tr]=Int_t(pptrk(1)); } @@ -931,34 +1113,29 @@ void AliMUON::Digitise(Int_t nev,Int_t bgr_ev,Option_t *option, Option_t *,Text_ } //end loop over clusters } // hit loop } // track loop - - //Int_t nentr1=list->GetEntriesFast(); - //printf(" \n counter, nentr1 %d %d\n",counter,nentr1); // open the file with background - if (Add ) { - ntracks =(Int_t)TrH1->GetEntries(); - //printf("background - icat,ntracks1 %d %d\n",icat,ntracks); - //printf("background - Start loop over tracks \n"); + if (addBackground) { + ntracks =(Int_t)fTrH1->GetEntries(); // // Loop over tracks // for (Int_t track=0; trackClear(); - if (fClusters2) fClusters2->Clear(); + if (fPadHits2) fPadHits2->Clear(); - TrH1->GetEvent(track); + fTrH1->GetEvent(track); // // Loop over hits - AliMUONhit* mHit; + AliMUONHit* mHit; for(int i=0;iGetEntriesFast();++i) - { - mHit=(AliMUONhit*) (*fHits2)[i]; + { + mHit=(AliMUONHit*) (*fHits2)[i]; Int_t nch = mHit->fChamber-1; // chamber number if (nch >9) continue; - iChamber = &(MUON->Chamber(nch)); + iChamber = &(pMUON->Chamber(nch)); Int_t rmin = (Int_t)iChamber->RInner(); Int_t rmax = (Int_t)iChamber->ROuter(); Float_t xbgr=mHit->fX; @@ -974,79 +1151,75 @@ void AliMUON::Digitise(Int_t nev,Int_t bgr_ev,Option_t *option, Option_t *,Text_ // // Loop over pad hits - for (AliMUONcluster* mPad= - (AliMUONcluster*)MUON->FirstPad(mHit,fClusters2); + for (AliMUONPadHit* mPad= + (AliMUONPadHit*)pMUON->FirstPad(mHit,fPadHits2); mPad; - mPad=(AliMUONcluster*)MUON->NextPad(fClusters2)) + mPad=(AliMUONPadHit*)pMUON->NextPad(fPadHits2)) { - + // mPad = (AliMUONPadHit*) (*fPadHits2)[j]; Int_t cathode = mPad->fCathode; // cathode number Int_t ipx = mPad->fPadX; // pad number on X Int_t ipy = mPad->fPadY; // pad number on Y - Int_t iqpad = Int_t(mPad->fQpad*kScale);// charge per pad -// Int_t iqpad = mPad->fQpad; // charge per pad + Int_t iqpad = Int_t(mPad->fQpad);// charge per pad if (cathode != (icat+1)) continue; - //if (!HitMap[nch]->CheckBoundary()) continue; - // fill the info array - Float_t thex, they; - segmentation=iChamber->GetSegmentationModel(cathode); - segmentation->GetPadCxy(ipx,ipy,thex,they); + Float_t thex, they, thez; + segmentation=iChamber->SegmentationModel(cathode); + segmentation->GetPadC(ipx,ipy,thex,they,thez); Float_t rpad=TMath::Sqrt(thex*thex+they*they); if (rpad < rmin || iqpad ==0 || rpad > rmax) continue; - - new((*p_adr)[countadr++]) TVector(2); - TVector &trinfo=*((TVector*) (*p_adr)[countadr-1]); - trinfo(0)=-1; // tag background - trinfo(1)=-1; - + new((*pAddress)[countadr++]) TVector(2); + TVector &trinfo=*((TVector*) (*pAddress)[countadr-1]); + trinfo(0)=-1; // tag background + trinfo(1)=-1; + digits[0]=ipx; digits[1]=ipy; digits[2]=iqpad; digits[3]=0; digits[4]=-1; - - AliMUONlist* pdigit; - // build the list of fired pads and update the info - if (!HitMap[nch]->TestHit(ipx, ipy)) { - list->AddAtAndExpand(new AliMUONlist(nch,digits),counter); - HitMap[nch]->SetHit(ipx, ipy, counter); + AliMUONTransientDigit* pdigit; + // build the list of fired pads and update the info + if (!hitMap[nch]->TestHit(ipx, ipy)) { + list->AddAtAndExpand(new AliMUONTransientDigit(nch,digits),counter); + + hitMap[nch]->SetHit(ipx, ipy, counter); counter++; - pdigit=(AliMUONlist*)list->At(list->GetLast()); + pdigit=(AliMUONTransientDigit*)list->At(list->GetLast()); // list of tracks - TObjArray *trlist=(TObjArray*)pdigit-> - TrackList(); - trlist->Add(&trinfo); + TObjArray *trlist=(TObjArray*)pdigit-> + TrackList(); + trlist->Add(&trinfo); } else { - pdigit=(AliMUONlist*) HitMap[nch]->GetHit(ipx, ipy); + pdigit=(AliMUONTransientDigit*) hitMap[nch]->GetHit(ipx, ipy); // update charge (*pdigit).fSignal+=iqpad; - + // update list of tracks - TObjArray* trlist=(TObjArray*)pdigit-> - TrackList(); - Int_t last_entry=trlist->GetLast(); - TVector *ptrk_p=(TVector*)trlist-> - At(last_entry); - TVector &ptrk=*ptrk_p; - Int_t last_track=Int_t(ptrk(0)); - if (last_track==-1) { - continue; - } else { - trlist->Add(&trinfo); - } + TObjArray* trlist=(TObjArray*)pdigit-> + TrackList(); + Int_t lastEntry=trlist->GetLast(); + TVector *pTrack=(TVector*)trlist-> + At(lastEntry); + TVector &ptrk=*pTrack; + Int_t lastTrack=Int_t(ptrk(0)); + if (lastTrack==-1) { + continue; + } else { + trlist->Add(&trinfo); + } // check the track list - Int_t nptracks=trlist->GetEntriesFast(); - if (nptracks > 0) { - for (Int_t tr=0;trAt(tr); - TVector &pptrk=*pptrk_p; - trk[tr]=Int_t(pptrk(0)); - chtrk[tr]=Int_t(pptrk(1)); - } - } // end if nptracks + Int_t nptracks=trlist->GetEntriesFast(); + if (nptracks > 0) { + for (Int_t tr=0;trAt(tr); + TVector &pptrk=*ppTrack; + trk[tr]=Int_t(pptrk(0)); + chtrk[tr]=Int_t(pptrk(1)); + } + } // end if nptracks } // end if pdigit } //end loop over clusters } // hit loop @@ -1054,114 +1227,100 @@ void AliMUON::Digitise(Int_t nev,Int_t bgr_ev,Option_t *option, Option_t *,Text_ //Int_t nentr2=list->GetEntriesFast(); //printf(" \n counter2, nentr2 %d %d \n",counter,nentr2); TTree *fAli=gAlice->TreeK(); - TFile *file; + TFile *file=NULL; if (fAli) file =fAli->GetCurrentFile(); file->cd(); - } // if Add + } // if addBackground + delete [] xhit; + delete [] yhit; Int_t tracks[10]; Int_t charges[10]; - //cout<<"start filling digits \n "<GetEntriesFast(); - //printf(" \n \n nentries %d \n",nentries); - // start filling the digits for (Int_t nent=0;nentAt(nent); + AliMUONTransientDigit *address=(AliMUONTransientDigit*)list->At(nent); if (address==0) continue; Int_t ich=address->fChamber; Int_t q=address->fSignal; - iChamber=(AliMUONchamber*) (*fChambers)[ich]; - AliMUONresponse * response=iChamber->GetResponseModel(); - Int_t adcmax= (Int_t) response->MaxAdc(); - // add white noise and do zero-suppression and signal truncation - Float_t MeanNoise = gRandom->Gaus(1, 0.2); - Float_t Noise = gRandom->Gaus(0, MeanNoise); - q+=(Int_t)Noise; - if (address->fPhysics !=0 ) address->fPhysics+=(Int_t)Noise; - if ( q <= zero_supm ) continue; - if ( q > adcmax) q=adcmax; + iChamber=(AliMUONChamber*) (*fChambers)[ich]; +// +// Digit Response (noise, threshold, saturation, ...) +// if (address->fPhysics !=0 ) address->fPhysics+=(Int_t)Noise; + AliMUONResponse * response=iChamber->ResponseModel(); + q=response->DigitResponse(q); + + if (!q) continue; + digits[0]=address->fPadX; digits[1]=address->fPadY; digits[2]=q; digits[3]=address->fPhysics; digits[4]=address->fHit; - //printf("fSignal, fPhysics fTrack %d %d %d \n",digits[2],digits[3],digits[4]); TObjArray* trlist=(TObjArray*)address->TrackList(); Int_t nptracks=trlist->GetEntriesFast(); //printf("nptracks, trlist %d %p\n",nptracks,trlist); - // this was changed to accomodate the real number of tracks - if (nptracks > 10) { - cout<<"Attention - nptracks > 10 "< 2) { - printf("Attention - nptracks > 2 %d \n",nptracks); - printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,digits[0],digits[1],q); - } - for (Int_t tr=0;trAt(tr); - if(!pp_p ) printf("pp_p - %p\n",pp_p); - TVector &pp =*pp_p; - tracks[tr]=Int_t(pp(0)); - charges[tr]=Int_t(pp(1)); + // this was changed to accomodate the real number of tracks + if (nptracks > 10) { + cout<<"Attention - nptracks > 10 "< 2) { + printf("Attention - nptracks > 2 %d \n",nptracks); + printf("cat,ich,ix,iy,q %d %d %d %d %d \n",icat,ich,digits[0],digits[1],q); + } + for (Int_t tr=0;trAt(tr); + if(!ppP ) printf("ppP - %p\n",ppP); + TVector &pp =*ppP; + tracks[tr]=Int_t(pp(0)); + charges[tr]=Int_t(pp(1)); //printf("tracks, charges - %d %d\n",tracks[tr],charges[tr]); - } //end loop over list of tracks for one pad + } //end loop over list of tracks for one pad // Sort list of tracks according to charge - if (nptracks > 1) { - SortTracks(tracks,charges,nptracks); - } - if (nptracks < 10 ) { - for (Int_t i=nptracks; i<10; i++) { - tracks[i]=0; - charges[i]=0; - } + if (nptracks > 1) { + SortTracks(tracks,charges,nptracks); + } + if (nptracks < 10 ) { + for (Int_t i=nptracks; i<10; i++) { + tracks[i]=0; + charges[i]=0; } - + } + // fill digits - MUON->AddDigits(ich,tracks,charges,digits); + pMUON->AddDigits(ich,tracks,charges,digits); + // delete trlist; } //cout<<"I'm out of the loops for digitisation"<GetEvent(nev); gAlice->TreeD()->Fill(); - TTree *TD=gAlice->TreeD(); - - Stat_t ndig=TD->GetEntries(); - cout<<"number of digits "<DigitsAddress(k); - int ndig=fDch->GetEntriesFast(); - printf (" i, ndig %d %d \n",k,ndig); - } - - MUON->ResetDigits(); + pMUON->ResetDigits(); list->Delete(); - for(Int_t ii=0;ii<10;++ii) { - if (HitMap[ii]) { - hm=HitMap[ii]; + + + for(Int_t ii=0;iiTreeD()->Write(hname); - // reset tree - gAlice->TreeD()->Reset(); - delete list; - //Int_t nadr=p_adr->GetEntriesFast(); - // printf(" \n \n nadr %d \n",nadr); - - p_adr->Clear(); - // gObjectTable->Print(); - + delete [] hitMap; + char hname[30]; + sprintf(hname,"TreeD%d",nev); + gAlice->TreeD()->Write(hname); + // reset tree + gAlice->TreeD()->Reset(); + delete list; + + pAddress->Delete(); + // gObjectTable->Print(); } void AliMUON::SortTracks(Int_t *tracks,Int_t *charges,Int_t ntr) @@ -1219,495 +1378,227 @@ void AliMUON::SortTracks(Int_t *tracks,Int_t *charges,Int_t ntr) } -void AliMUON::FindClusters(Int_t nev,Int_t last_entry) -{ - -// -// Loop on chambers and on cathode planes -// - for (Int_t icat=0;icat<2;icat++) { - gAlice->ResetDigits(); - gAlice->TreeD()->GetEvent(last_entry+icat); // spurious +1 ... - if (nev < 10) printf("last_entry , icat - %d %d \n",last_entry,icat); - //gAlice->TreeD()->GetEvent(icat+1); // spurious +1 ... - - for (Int_t ich=0;ich<10;ich++) { - AliMUONchamber* iChamber=(AliMUONchamber*) (*fChambers)[ich]; - TClonesArray *MUONdigits = this->DigitsAddress(ich); - if (MUONdigits == 0) continue; - // - // Get ready the current chamber stuff - // - AliMUONresponse* response = iChamber->GetResponseModel(); - AliMUONsegmentation* seg = iChamber->GetSegmentationModel(icat+1); - AliMUONClusterFinder* rec = iChamber->GetReconstructionModel(); - //printf("icat, ich, seg - %d %d %p\n",icat,ich,seg); - if (seg) { - rec->SetSegmentation(seg); - rec->SetResponse(response); - rec->SetDigits(MUONdigits); - rec->SetChamber(ich); - if (nev==0) rec->CalibrateCOG(); - rec->FindRawClusters(); - } - //printf("Finish FindRawClusters for cathode %d in chamber %d\n",icat,ich); - - TClonesArray *fRch; - fRch=RawClustAddress(ich); - fRch->Sort(); - // it seems to work - - - } // for ich - // fill the tree - TTree *TR=gAlice->TreeR(); - - gAlice->TreeR()->Fill(); - - Stat_t nent=TR->GetEntries(); - cout<<"number of entries "<GetEntriesFast(); - printf (" i, nraw %d %d \n",i,nraw); +//___________________________________________ +void AliMUON::Trigger(Int_t nev){ +// call the Trigger Algorithm and fill TreeR + + Int_t singlePlus[3] = {0,0,0}; + Int_t singleMinus[3] = {0,0,0}; + Int_t singleUndef[3] = {0,0,0}; + Int_t pairUnlike[3] = {0,0,0}; + Int_t pairLike[3] = {0,0,0}; + + ResetTrigger(); + AliMUONTriggerDecision* decision= new AliMUONTriggerDecision(1); + decision->Trigger(); + decision->GetGlobalTrigger(singlePlus, singleMinus, singleUndef, + pairUnlike, pairLike); +// add a local trigger in the list + AddGlobalTrigger(singlePlus, singleMinus, singleUndef, pairUnlike, pairLike); + Int_t i; + + for (Int_t icirc=0; icircGetITrigger(icirc)==1) { + Int_t localtr[7]={0,0,0,0,0,0,0}; + Int_t loLpt[2]={0,0}; Int_t loHpt[2]={0,0}; Int_t loApt[2]={0,0}; + decision->GetLutOutput(icirc, loLpt, loHpt, loApt); + localtr[0] = icirc; + localtr[1] = decision->GetStripX11(icirc); + localtr[2] = decision->GetDev(icirc); + localtr[3] = decision->GetStripY11(icirc); + for (i=0; i<2; i++) { // convert the Lut output in 1 digit + localtr[4] = localtr[4]+Int_t(loLpt[i]*TMath::Power(2,i)); + localtr[5] = localtr[5]+Int_t(loHpt[i]*TMath::Power(2,i)); + localtr[6] = localtr[6]+Int_t(loApt[i]*TMath::Power(2,i)); + } + AddLocalTrigger(localtr); // add a local trigger in the list } - ResetRawClusters(); - - } // for icat + } + delete decision; + gAlice->TreeR()->Fill(); + ResetTrigger(); char hname[30]; sprintf(hname,"TreeR%d",nev); gAlice->TreeR()->Write(hname); gAlice->TreeR()->Reset(); - - //gObjectTable->Print(); - + printf("\n End of trigger for event %d", nev); } - -//______________________________________________________________________________ -//_____________________________________________________________________________ -void AliMUON::CathodeCorrelation(Int_t nev) -{ -// Correlates the clusters on the two cathode planes and build a list of -// other possible combinations (potential ghosts) - for the moment use the -// criteria of minimum distance between the CoGs of the two correlated -// clusters - -// -// Loop on chambers and on clusters on the cathode plane with the highest -// number of clusters - - static Bool_t first=kTRUE; - - AliMUONRawCluster *mRaw1; - AliMUONRawCluster *mRaw2; - AliMUONchamber *iChamber; - AliMUONsegmentation *seg; - TArrayF x1, y1, x2, y2, q1, q2; - x1.Set(5000); - x2.Set(5000); - y1.Set(5000); - y2.Set(5000); - q1.Set(5000); - q2.Set(5000); - -// Get pointers to Alice detectors and Digits containers - TTree *TR = gAlice->TreeR(); - Int_t nent=(Int_t)TR->GetEntries(); - if (nev < 10) printf("Found %d entries in the tree (must be one per cathode per event! + 1empty)\n",nent); - - - Int_t idx[4]; - Float_t xc2[4],yc2[4]; - Float_t xrec2, yrec2; - Float_t xd0, xdif, ydif; - Float_t ysrch,xd,xmax,ymax; - Int_t ilow, iup, iraw1, i; - // - Float_t xarray[50]; - Float_t xdarray[50]; - Float_t yarray[50]; - Float_t qarray[50]; - Int_t idx2[50]; - - // Int_t nraw[2], entry,cathode; - - for (i=0;i<50;i++) { - xdarray[i]=1100.; - xarray[i]=0.; - yarray[i]=0.; - qarray[i]=0.; - idx2[i]=-1; - } - for (i=0;i<4;i++) { - idx[i]=-1; - xc2[i]=0.; - yc2[i]=0.; - } - - // access to the Raw Clusters tree - for (Int_t ich=0;ich<10;ich++) { - iChamber = &(Chamber(ich)); - TClonesArray *MUONrawclust = RawClustAddress(ich); - ResetRawClusters(); - TR->GetEvent(nent-2); - //TR->GetEvent(1); - Int_t nrawcl1 = MUONrawclust->GetEntries(); - // printf("Found %d raw clusters for cathode 1 in chamber %d \n" - // ,nrawcl1,ich+1); - if (!nrawcl1) continue; - - seg = iChamber->GetSegmentationModel(1); - // loop over raw clusters of first cathode - for (iraw1=0; iraw1UncheckedAt(iraw1); - x1[iraw1]=mRaw1->fX; - y1[iraw1]=mRaw1->fY; - q1[iraw1]=(Float_t)mRaw1->fQ; //maybe better fPeakSignal - } // rawclusters cathode 1 +//____________________________________________ +void AliMUON::FindClusters(Int_t nev,Int_t lastEntry) +{ + TClonesArray *dig1, *dig2; + Int_t ndig, k; + dig1 = new TClonesArray("AliMUONDigit",1000); + dig2 = new TClonesArray("AliMUONDigit",1000); + AliMUONDigit *digit; // - // Get information from 2nd cathode - ResetRawClusters(); - TR->GetEvent(nent-1); - //TR->GetEvent(2); - Int_t nrawcl2 = MUONrawclust->GetEntries(); - if (!nrawcl2) { - for (iraw1=0; iraw1UncheckedAt(iraw2); - x2[iraw2]=mRaw2->fX; - y2[iraw2]=mRaw2->fY; - q2[iraw2]=(Float_t)mRaw2->fQ; - } // rawclusters cathode 2 +// Loop on chambers and on cathode planes // -// Initalisation finished - for (iraw1=0; iraw1GetPadIxy(x1[iraw1],y1[iraw1],ix,iy); - Int_t isec=seg->Sector(ix,iy); - // range to look for ghosts ?! - if (ich < 5) { - ymax = seg->Dpy(isec)*7/2; - xmax = seg->Dpx(isec)*7/2; - } else { - ymax = seg->Dpy(isec)*13/2; - xmax = seg->Dpx(isec)*3/2; - } - ysrch=ymax+y1[iraw1]; - - ilow = AliMUONRawCluster:: - BinarySearch(ysrch-2*ymax,y2,0,nrawcl2+1); - iup= AliMUONRawCluster:: - BinarySearch(ysrch,y2,ilow,nrawcl2+1); - if (ilow<0 || iup <0 || iup>nrawcl2) continue; - Int_t counter=0; - for (Int_t iraw2=ilow; iraw2<=iup; iraw2++) { - xrec2=x2[iraw2]; - yrec2=y2[iraw2]; - xdif=x1[iraw1]-xrec2; - ydif=y1[iraw1]-yrec2; - xd=TMath::Sqrt(xdif*xdif+ydif*ydif); - if (iraw2==ilow) { - if (ilow==iup) - xd0=TMath:: - Sqrt(2*xmax*2*xmax+2*ymax*2*ymax); - else xd0=101.; - } - Float_t qdif=TMath::Abs(q1[iraw1]-q2[iraw2])/q1[iraw1]; - - if (x1[iraw1]*xrec2 > 0) { - if (xd <= xd0 ) { -// printf("q1, q2 qdif % f %f %f \n",q1[iraw1],q2[iraw2],qdif); -// printf("x1, x2 y1 y2 % f %f %f %f \n",x1[iraw1],xrec2,y1[iraw1],yrec2); - //if (qdif <0.3) { //check this number - - xd0=xd; - idx2[counter]=iraw2; - xdarray[counter]=xd; - xarray[counter]=xdif; - yarray[counter]=ydif; - qarray[counter]=qdif; - counter++; - // } - - } - } // check for same quadrant - } // loop over 2nd cathode range - - - if (counter >=2) { - AliMUONRawCluster:: - SortMin(idx2,xdarray,xarray,yarray,qarray,counter); - if (xdarray[0]Dpx(isec) && xdarray[1]Dpx(isec)) { - if (qarray[0]>qarray[1]){ - Int_t swap=idx2[0]; - idx2[0]=idx2[1]; - idx2[1]=swap; - } - } - } - int imax; - if (counter <3) imax=counter; - else imax=3; - - for (int i=0;i= 0 && idx2[i] < nrawcl2) { - if (xarray[i] > xmax || yarray[i] > 2*ymax) - continue; - idx[i]=idx2[i]; - xc2[i]=x2[idx2[i]]; - yc2[i]=y2[idx2[i]]; - } - } - // add info about the cluster on the 'starting' cathode - - idx[3]=iraw1; - xc2[3]=x1[iraw1]; - yc2[3]=y1[iraw1]; - //if (idx[0] <0) printf("iraw1 imax idx2[0] idx[0] %d %d %d %d\n",iraw1,imax,idx2[0],idx[0]); - AddCathCorrel(ich,idx,xc2,yc2); - // reset - for (Int_t ii=0;iiFill(); - //Int_t nentries=(Int_t)TC->GetEntries(); - //cout<<"number entries in tree of correlated clusters "<GetEntriesFast(); - printf (" ii, ncor %d %d \n",ii,ncor); - if (ncor>=2) countch++; - } - - // write - char hname[30]; - sprintf(hname,"TreeC%d",nev); - TC->Write(hname); - // reset tree - ResetCorrelation(); - TC->Reset(); - - if (countch==10) countev++; - printf("countev - %d\n",countev); -// gObjectTable->Print(); - - -} - - -//_____________________________________________________________________________ - -void AliMUON::MakeTreeC(Option_t *option) -{ - char *C = strstr(option,"C"); - if (C && !fTreeC) fTreeC = new TTree("TC","CathodeCorrelation"); - -// Create a branch for correlation - - const Int_t buffersize = 4000; - char branchname[30]; - -// one branch for correlation per chamber - for (int i=0; i<10 ;i++) { - sprintf(branchname,"%sCorrelation%d",GetName(),i+1); - - if (fCathCorrel && fTreeC) { - TreeC()->Branch(branchname,&((*fCathCorrel)[i]), buffersize); - printf("Making Branch %s for correlation in chamber %d\n",branchname,i+1); - } - } -} - -//_____________________________________________________________________________ -void AliMUON::GetTreeC(Int_t event) -{ - - // set the branch address - char treeName[20]; - char branchname[30]; - - ResetCorrelation(); - if (fTreeC) { - delete fTreeC; - } - - sprintf(treeName,"TreeC%d",event); - fTreeC = (TTree*)gDirectory->Get(treeName); - - - TBranch *branch; - if (fTreeC) { - for (int i=0; i<10; i++) { - sprintf(branchname,"%sCorrelation%d",GetName(),i+1); - if (fCathCorrel) { - branch = fTreeC->GetBranch(branchname); - if (branch) branch->SetAddress(&((*fCathCorrel)[i])); - } + for (Int_t ich=0;ich<10;ich++) { + AliMUONChamber* iChamber=(AliMUONChamber*) (*fChambers)[ich]; + AliMUONClusterFinderVS* rec = iChamber->ReconstructionModel(); + gAlice->ResetDigits(); + gAlice->TreeD()->GetEvent(lastEntry); + TClonesArray *muonDigits = this->DigitsAddress(ich); + ndig=muonDigits->GetEntriesFast(); + printf("\n 1 Found %d digits in %p %d", ndig, muonDigits,ich); + TClonesArray &lhits1 = *dig1; + Int_t n=0; + for (k=0; kUncheckedAt(k); + if (rec->TestTrack(digit->fTracks[0])) + new(lhits1[n++]) AliMUONDigit(*digit); + } + gAlice->ResetDigits(); + gAlice->TreeD()->GetEvent(lastEntry+1); + muonDigits = this->DigitsAddress(ich); + ndig=muonDigits->GetEntriesFast(); + printf("\n 2 Found %d digits in %p %d", ndig, muonDigits, ich); + TClonesArray &lhits2 = *dig2; + n=0; + + for (k=0; kUncheckedAt(k); + if (rec->TestTrack(digit->fTracks[0])) + new(lhits2[n++]) AliMUONDigit(*digit); } - } else { - printf("ERROR: cannot find CathodeCorrelation Tree for event:%d\n",event); - } - - // gObjectTable->Print(); + if (rec) { + AliMUONClusterInput::Instance()->SetDigits(ich, dig1, dig2); + rec->FindRawClusters(); + } + dig1->Delete(); + dig2->Delete(); + } // for ich + gAlice->TreeR()->Fill(); + ResetRawClusters(); + char hname[30]; + sprintf(hname,"TreeR%d",nev); + gAlice->TreeR()->Write(hname); + gAlice->TreeR()->Reset(); + printf("\n End of cluster finding for event %d", nev); + + delete dig1; + delete dig2; + //gObjectTable->Print(); } - + void AliMUON::Streamer(TBuffer &R__b) { // Stream an object of class AliMUON. - AliMUONchamber *iChamber; - AliMUONsegmentation *segmentation; - AliMUONresponse *response; - TClonesArray *digitsaddress; - TClonesArray *rawcladdress; - TClonesArray *corcladdress; - // TObjArray *clustaddress; - - if (R__b.IsReading()) { - Version_t R__v = R__b.ReadVersion(); if (R__v) { } - AliDetector::Streamer(R__b); - R__b >> fNclusters; - R__b >> fClusters; // diff - R__b >> fDchambers; - R__b >> fRawClusters; - R__b >> fCathCorrel; - R__b.ReadArray(fNdch); - R__b.ReadArray(fNrawch); - R__b.ReadArray(fNcorch); - // - R__b >> fAccCut; - R__b >> fAccMin; - R__b >> fAccMax; - // - // modifs perso - R__b >> fSPxzCut; - R__b >> fSSigmaCut; - R__b >> fSXPrec; - R__b >> fSYPrec; - // - R__b >> fChambers; + AliMUONChamber *iChamber; + AliMUONTriggerCircuit *iTriggerCircuit; + AliSegmentation *segmentation; + AliMUONResponse *response; + TClonesArray *digitsaddress; + TClonesArray *rawcladdress; + Int_t i; + if (R__b.IsReading()) { + Version_t R__v = R__b.ReadVersion(); if (R__v) { } + AliDetector::Streamer(R__b); + R__b >> fNPadHits; + R__b >> fPadHits; // diff + R__b >> fNLocalTrigger; + R__b >> fLocalTrigger; + R__b >> fNGlobalTrigger; + R__b >> fGlobalTrigger; + R__b >> fDchambers; + R__b >> fRawClusters; + R__b.ReadArray(fNdch); + R__b.ReadArray(fNrawch); + R__b >> fAccCut; + R__b >> fAccMin; + R__b >> fAccMax; + R__b >> fChambers; + R__b >> fTriggerCircuits; + for (i =0; iStreamer(R__b); + } // Stream chamber related information - for (Int_t i =0; i<10; i++) { - iChamber=(AliMUONchamber*) (*fChambers)[i]; - iChamber->Streamer(R__b); - if (iChamber->Nsec()==1) { - segmentation=iChamber->GetSegmentationModel(1); - segmentation->Streamer(R__b); - } else { - segmentation=iChamber->GetSegmentationModel(1); - segmentation->Streamer(R__b); - segmentation=iChamber->GetSegmentationModel(2); - segmentation->Streamer(R__b); + for (i =0; iStreamer(R__b); + if (iChamber->Nsec()==1) { + segmentation=iChamber->SegmentationModel(1); + if (segmentation) + segmentation->Streamer(R__b); + } else { + segmentation=iChamber->SegmentationModel(1); + if (segmentation) + segmentation->Streamer(R__b); + if (segmentation) + segmentation=iChamber->SegmentationModel(2); + segmentation->Streamer(R__b); + } + response=iChamber->ResponseModel(); + if (response) + response->Streamer(R__b); + digitsaddress=(TClonesArray*) (*fDchambers)[i]; + digitsaddress->Streamer(R__b); + if (i < AliMUONConstants::NTrackingCh()) { + rawcladdress=(TClonesArray*) (*fRawClusters)[i]; + rawcladdress->Streamer(R__b); + } } - response=iChamber->GetResponseModel(); - response->Streamer(R__b); - digitsaddress=(TClonesArray*) (*fDchambers)[i]; - digitsaddress->Streamer(R__b); - rawcladdress=(TClonesArray*) (*fRawClusters)[i]; - rawcladdress->Streamer(R__b); - corcladdress=(TClonesArray*) (*fCathCorrel)[i]; - corcladdress->Streamer(R__b); - } - - } else { - R__b.WriteVersion(AliMUON::IsA()); - AliDetector::Streamer(R__b); - R__b << fNclusters; - R__b << fClusters; // diff - R__b << fDchambers; - R__b << fRawClusters; - R__b << fCathCorrel; - R__b.WriteArray(fNdch, 10); - R__b.WriteArray(fNrawch, 10); - R__b.WriteArray(fNcorch, 10); - // - R__b << fAccCut; - R__b << fAccMin; - R__b << fAccMax; - // - // modifs perso - R__b << fSPxzCut; - R__b << fSSigmaCut; - R__b << fSXPrec; - R__b << fSYPrec; - // - R__b << fChambers; -// Stream chamber related information - for (Int_t i =0; i<10; i++) { - iChamber=(AliMUONchamber*) (*fChambers)[i]; - iChamber->Streamer(R__b); - if (iChamber->Nsec()==1) { - segmentation=iChamber->GetSegmentationModel(1); - segmentation->Streamer(R__b); - } else { - segmentation=iChamber->GetSegmentationModel(1); - segmentation->Streamer(R__b); - segmentation=iChamber->GetSegmentationModel(2); - segmentation->Streamer(R__b); + + } else { + R__b.WriteVersion(AliMUON::IsA()); + AliDetector::Streamer(R__b); + R__b << fNPadHits; + R__b << fPadHits; // diff + R__b << fNLocalTrigger; + R__b << fLocalTrigger; + R__b << fNGlobalTrigger; + R__b << fGlobalTrigger; + R__b << fDchambers; + R__b << fRawClusters; + R__b.WriteArray(fNdch, AliMUONConstants::NCh()); + R__b.WriteArray(fNrawch, AliMUONConstants::NTrackingCh()); + + R__b << fAccCut; + R__b << fAccMin; + R__b << fAccMax; + + R__b << fChambers; + R__b << fTriggerCircuits; + for (i =0; iStreamer(R__b); + } + for (i =0; iStreamer(R__b); + if (iChamber->Nsec()==1) { + segmentation=iChamber->SegmentationModel(1); + if (segmentation) + segmentation->Streamer(R__b); + } else { + segmentation=iChamber->SegmentationModel(1); + if (segmentation) + segmentation->Streamer(R__b); + segmentation=iChamber->SegmentationModel(2); + if (segmentation) + segmentation->Streamer(R__b); + } + response=iChamber->ResponseModel(); + if (response) + response->Streamer(R__b); + digitsaddress=(TClonesArray*) (*fDchambers)[i]; + digitsaddress->Streamer(R__b); + if (i < AliMUONConstants::NTrackingCh()) { + rawcladdress=(TClonesArray*) (*fRawClusters)[i]; + rawcladdress->Streamer(R__b); + } } - response=iChamber->GetResponseModel(); - response->Streamer(R__b); - digitsaddress=(TClonesArray*) (*fDchambers)[i]; - digitsaddress->Streamer(R__b); - rawcladdress=(TClonesArray*) (*fRawClusters)[i]; - rawcladdress->Streamer(R__b); - corcladdress=(TClonesArray*) (*fCathCorrel)[i]; - corcladdress->Streamer(R__b); } - } } -AliMUONcluster* AliMUON::FirstPad(AliMUONhit* hit, TClonesArray *clusters) +AliMUONPadHit* AliMUON::FirstPad(AliMUONHit* hit, TClonesArray *clusters) { // // Initialise the pad iterator @@ -1715,1493 +1606,47 @@ AliMUONcluster* AliMUON::FirstPad(AliMUONhit* hit, TClonesArray *clusters) TClonesArray *theClusters = clusters; Int_t nclust = theClusters->GetEntriesFast(); if (nclust && hit->fPHlast > 0) { - sMaxIterPad=hit->fPHlast; - sCurIterPad=hit->fPHfirst; - return (AliMUONcluster*) clusters->UncheckedAt(sCurIterPad-1); + AliMUON::fMaxIterPad=hit->fPHlast; + AliMUON::fCurIterPad=hit->fPHfirst; + return (AliMUONPadHit*) clusters->UncheckedAt(AliMUON::fCurIterPad-1); } else { return 0; } } -AliMUONcluster* AliMUON::NextPad(TClonesArray *clusters) +AliMUONPadHit* AliMUON::NextPad(TClonesArray *clusters) { - sCurIterPad++; - if (sCurIterPad <= sMaxIterPad) { - return (AliMUONcluster*) clusters->UncheckedAt(sCurIterPad-1); + AliMUON::fCurIterPad++; + if (AliMUON::fCurIterPad <= AliMUON::fMaxIterPad) { + return (AliMUONPadHit*) clusters->UncheckedAt(AliMUON::fCurIterPad-1); } else { return 0; } } -//////////////////////////// modifs perso /////////////// - -static TTree *ntuple_global; -static TFile *hfile_global; - -// variables of the tracking ntuple -struct { - Int_t ievr; // number of event - Int_t ntrackr; // number of tracks per event - Int_t istatr[500]; // 1 = good muon, 2 = ghost, 0 = something else - Int_t isignr[500]; // sign of the track - Float_t pxr[500]; // x momentum of the reconstructed track - Float_t pyr[500]; // y momentum of the reconstructed track - Float_t pzr[500]; // z momentum of the reconstructed track - Float_t zvr[500]; // z vertex - Float_t chi2r[500]; // chi2 of the fit of the track with the field map - Float_t pxv[500]; // x momentum at vertex - Float_t pyv[500]; // y momentum at vertex - Float_t pzv[500]; // z momentum at vertex -} ntuple_st; AliMUONRawCluster *AliMUON::RawCluster(Int_t ichamber, Int_t icathod, Int_t icluster) { - TClonesArray *MUONrawclust = RawClustAddress(ichamber); + TClonesArray *muonRawCluster = RawClustAddress(ichamber); ResetRawClusters(); - TTree *TR = gAlice->TreeR(); - Int_t nent=(Int_t)TR->GetEntries(); - TR->GetEvent(nent-2+icathod-1); - //TR->GetEvent(icathod); - //Int_t nrawcl = (Int_t)MUONrawclust->GetEntriesFast(); + TTree *treeR = gAlice->TreeR(); + Int_t nent=(Int_t)treeR->GetEntries(); + treeR->GetEvent(nent-2+icathod-1); + //treeR->GetEvent(icathod); + //Int_t nrawcl = (Int_t)muonRawCluster->GetEntriesFast(); - AliMUONRawCluster * mRaw = (AliMUONRawCluster*)MUONrawclust->UncheckedAt(icluster); + AliMUONRawCluster * mRaw = (AliMUONRawCluster*)muonRawCluster->UncheckedAt(icluster); //printf("RawCluster _ nent nrawcl icluster mRaw %d %d %d%p\n",nent,nrawcl,icluster,mRaw); return mRaw; } -void AliMUON::Reconst(Int_t &ifit, Int_t &idebug, Int_t bgd_ev, Int_t &nev, Int_t &idres, Int_t &ireadgeant, Option_t *option,Text_t *filename) -{ - // - // open kine and hits tree of background file for reconstruction of geant hits - // call tracking fortran program - static Bool_t first=kTRUE; - static TFile *File; - char *Add = strstr(option,"Add"); - - if (Add ) { // only in case of background with geant hits - if(first) { - fFileName=filename; - cout<<"filename "<cd(); - if(fHits2) fHits2->Clear(); - if(fParticles2) fParticles2->Clear(); - if(TrH1) delete TrH1; - TrH1=0; - if(TK1) delete TK1; - TK1=0; - // Get Hits Tree header from file - char treeName[20]; - sprintf(treeName,"TreeH%d",bgd_ev); - TrH1 = (TTree*)gDirectory->Get(treeName); - if (!TrH1) { - printf("ERROR: cannot find Hits Tree for event:%d\n",bgd_ev); - } - // set branch addresses - TBranch *branch; - char branchname[30]; - sprintf(branchname,"%s",GetName()); - if (TrH1 && fHits2) { - branch = TrH1->GetBranch(branchname); - if (branch) branch->SetAddress(&fHits2); - } - TrH1->GetEntries(); - // get the Kine tree - sprintf(treeName,"TreeK%d",bgd_ev); - TK1 = (TTree*)gDirectory->Get(treeName); - if (!TK1) { - printf("ERROR: cannot find Kine Tree for event:%d\n",bgd_ev); - } - // set branch addresses - if (TK1) - TK1->SetBranchAddress("Particles", &fParticles2); - TK1->GetEvent(0); - - // get back to the first file - TTree *TK = gAlice->TreeK(); - TFile *file1 = 0; - if (TK) file1 = TK->GetCurrentFile(); - file1->cd(); - - } // end if Add - - // call tracking fortran program - reconstmuon(ifit,idebug,nev,idres,ireadgeant); -} - - -void AliMUON::InitTracking(Double_t &seff, Double_t &sb0, Double_t &sbl3) -{ - // - // introduce in fortran program somme parameters and cuts for tracking - // create output file "reconst.root" (histos + ntuple) - cutpxz(fSPxzCut); // Pxz cut (GeV/c) to begin the track finding - sigmacut(fSSigmaCut); // Number of sigmas delimiting the searching areas - xpreci(fSXPrec); // Chamber precision in X (cm) - ypreci(fSYPrec); // Chamber precision in Y (cm) - reco_init(seff,sb0,sbl3); -} - -void AliMUON::FinishEvent() -{ - TTree *TK = gAlice->TreeK(); - if (TK) { - TFile *file1 = TK->GetCurrentFile(); - if(file1) file1->cd(); - } -} - -void AliMUON::CloseTracking() -{ - // - // write histos and ntuple to "reconst.root" file - reco_term(); -} - -void chfill(Int_t &id, Float_t &x, Float_t &, Float_t &) -{ - // - // fill histo like hfill in fortran - char name[5]; - sprintf(name,"h%d",id); - TH1F *h1 = (TH1F*) gDirectory->Get(name); - h1->Fill(x); -} - -void chfill2(Int_t &id, Float_t &x, Float_t &y, Float_t &w) -{ - // - // fill histo like hfill2 in fortran - char name[5]; - sprintf(name,"h%d",id); - TH2F *h2 = (TH2F*) gDirectory->Get(name); - h2->Fill(x,y,w); -} - -void chf1(Int_t &id, Float_t &x, Float_t &w) -{ - // - // fill histo like hf1 in fortran - char name[5]; - sprintf(name,"h%d",id); - TH1F *h1 = (TH1F*) gDirectory->Get(name); - h1->Fill(x,w); -} - -void hist_create() -{ - // - // Create an output file ("reconst.root") - // Create some histograms and an ntuple - - hfile_global = new TFile("reconst.root","RECREATE","Ntuple - reconstruction"); - - ntuple_global = new TTree("ntuple","Reconst ntuple"); - ntuple_global->Branch("ievr",&ntuple_st.ievr,"ievr/I"); - ntuple_global->Branch("ntrackr",&ntuple_st.ntrackr,"ntrackr/I"); - ntuple_global->Branch("istatr",&ntuple_st.istatr[0],"istatr[500]/I"); - ntuple_global->Branch("isignr",&ntuple_st.isignr[0],"isignr[500]/I"); - ntuple_global->Branch("pxr",&ntuple_st.pxr[0],"pxr[500]/F"); - ntuple_global->Branch("pyr",&ntuple_st.pyr[0],"pyr[500]/F"); - ntuple_global->Branch("pzr",&ntuple_st.pzr[0],"pzr[500]/F"); - ntuple_global->Branch("zvr",&ntuple_st.zvr[0],"zvr[500]/F"); - ntuple_global->Branch("chi2r",&ntuple_st.chi2r[0],"chi2r[500]/F"); - ntuple_global->Branch("pxv",&ntuple_st.pxv[0],"pxv[500]/F"); - ntuple_global->Branch("pyv",&ntuple_st.pyv[0],"pyv[500]/F"); - ntuple_global->Branch("pzv",&ntuple_st.pzv[0],"pzv[500]/F"); - - // test aliroot - - new TH1F("h100","particule id du hit geant",20,0.,20.); - new TH1F("h101","position en x du hit geant",100,-200.,200.); - new TH1F("h102","position en y du hit geant",100,-200.,200.); - new TH1F("h103","chambre de tracking concernee",15,0.,14.); - new TH1F("h104","moment ptot du hit geant",50,0.,100.); - new TH1F("h105","px au vertex",50,0.,20.); - new TH1F("h106","py au vertex",50,0.,20.); - new TH1F("h107","pz au vertex",50,0.,20.); - new TH1F("h108","position zv",50,-15.,15.); - new TH1F("h109","position en x du hit reconstruit",100,-300.,300.); - new TH1F("h110","position en y du hit reconstruit",100,-300.,300.); - new TH1F("h111","delta x ",100,-0.4,0.4); - new TH1F("h112","delta y ",100,-0.4,0.4); - - char hname[30]; - char hname1[30]; - for (int i=0;i<10;i++) { - sprintf(hname,"deltax%d",i); - sprintf(hname1,"h12%d",i); - new TH1F(hname1,hname ,100,-0.4,0.4); - sprintf(hname,"deltay%d",i); - sprintf(hname1,"h13%d",i); - new TH1F(hname1,hname ,100,-0.4,0.4); - } - new TH2F("h2000","VAR X st. 5",30,3.0,183.0,100,0.,25.); - new TH2F("h2001","VAR Y st. 5",30,3.0,183.0,100,0.,25.); - - new TH2F("h2500","P vs X HHIT",30,3.0,183.0,200,0.,200.); - new TH2F("h2501","P vs X HHIT**2",30,3.0,183.0,200,0.,5000.); - new TH2F("h2502","P vs X EPH2 st. 5",30,3.0,183.0,100,0.,0.000005); - new TH2F("h2503","P vs X EAL2 st. 5",30,3.0,183.0,100,0.,0.01); - //new TH2F("h2504","P vs X EXM2 st. 5",30,3.0,183.0,100,0.,1.5); - new TH2F("h2504","P vs X EXM2 st. 5",30,3.0,183.0,100,0.,0.1); - new TH2F("h2505","P vs X EYM2 st. 5",30,3.0,183.0,100,0.,30.); - - new TH2F("h2507","P vs X EPH st. 5",30,3.0,183.0,100,0.,0.003); - new TH2F("h2508","P vs X EAL st. 5",30,3.0,183.0,100,0.,0.3); - //new TH2F("h2509","P vs X EXM st. 5",30,3.0,183.0,100,0.,1.5); - new TH2F("h2509","P vs X EXM st. 5",30,3.0,183.0,100,0.,0.4); - new TH2F("h2510","P vs X EYM st. 5",30,3.0,183.0,100,0.,30.); - - new TH2F("h2511","P vs X EPH cut st. 5",30,3.0,183.0,100,0.,0.01); - new TH2F("h2512","P vs X EAL cut st. 5",30,3.0,183.0,100,0.,0.3); - //new TH2F("h2513","P vs X EXM cut st. 5",30,3.0,183.0,100,0.,1.5); - new TH2F("h2513","P vs X EXM cut st. 5",30,3.0,183.0,100,0.,0.4); - new TH2F("h2514","P vs X EYM cut st. 5",30,3.0,183.0,100,0.,30.); - // 4 - new TH2F("h2400","P vs X HHIT",30,3.0,183.0,200,0.,200.); - new TH2F("h2401","P vs X HHIT**2",30,3.0,183.0,200,0.,5000.); - new TH2F("h2402","P vs X EPH2 st. 4",30,3.0,183.0,100,0.,0.000005); - new TH2F("h2403","P vs X EAL2 st. 4",30,3.0,183.0,100,0.,0.05); - //new TH2F("h2404","P vs X EXM2 st. 4",30,3.0,183.0,100,0.,1.5); - new TH2F("h2404","P vs X EXM2 st. 4",30,3.0,183.0,100,0.,0.1); - new TH2F("h2405","P vs X EYM2 st. 4",30,3.0,183.0,100,0.,30.); - - new TH2F("h2407","P vs X EPH st. 4",30,3.0,183.0,100,0.,0.003); - new TH2F("h2408","P vs X EAL st. 4",30,3.0,183.0,100,0.,0.3); - //new TH2F("h2409","P vs X EXM st. 4",30,3.0,183.0,100,0.,1.5); - new TH2F("h2409","P vs X EXM st. 4",30,3.0,183.0,100,0.,0.1); - new TH2F("h2410","P vs X EYM st. 4",30,3.0,183.0,100,0.,30.); - - new TH2F("h2411","P vs X EPH cut st. 4",30,3.0,183.0,100,0.,0.01); - new TH2F("h2412","P vs X EAL cut st. 4",30,3.0,183.0,100,0.,0.3); - //new TH2F("h2413","P vs X EXM cut st. 4",30,3.0,183.0,100,0.,1.5); - new TH2F("h2413","P vs X EXM cut st. 4",30,3.0,183.0,100,0.,0.1); - new TH2F("h2414","P vs X EYM cut st. 4",30,3.0,183.0,100,0.,30.); - // 3 - new TH1F("h2301","P2",30,3.0,183.0); - new TH2F("h2302","P2 vs X EPH2 st. 3",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2303","P2 vs X EAL2 st. 3",30,3.0,183.0,100,0.,0.0005); - //new TH2F("h2304","P2 vs X EXM2 st. 3",30,3.0,183.0,100,0.,1.5); - new TH2F("h2304","P2 vs X EXM2 st. 3",30,3.0,183.0,100,0.,2.); - new TH2F("h2305","P2 vs X EYM2 st. 3",30,3.0,183.0,100,0.,3.); - - new TH2F("h2307","P vs X EPH2 st. 3",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2308","P vs X EAL2 st. 3",30,3.0,183.0,100,0.,0.005); - //new TH2F("h2309","P vs X EXM2 st. 3",30,3.0,183.0,100,0.,1.5); - new TH2F("h2309","P vs X EXM2 st. 3",30,3.0,183.0,100,0.,2.); - new TH2F("h2310","P vs X EYM2 st. 3",30,3.0,183.0,100,0.,3.); - - new TH2F("h2311","P vs X EPH cut st. 3",30,3.0,183.0,100,0.,0.06); - new TH2F("h2312","P vs X EAL cut st. 3",30,3.0,183.0,100,0.,0.05); - //new TH2F("h2313","P vs X EXM cut st. 3",30,3.0,183.0,100,0.,1.5); - new TH2F("h2313","P vs X EXM cut st. 3",30,3.0,183.0,100,0.,6.); - new TH2F("h2314","P vs X EYM cut st. 3",30,3.0,183.0,100,0.,7.); - - new TH2F("h2315","P2 vs X EPH cut st. 3",30,3.0,183.0,100,0.,0.06); - new TH2F("h2316","P2 vs X EAL cut st. 3",30,3.0,183.0,100,0.,0.05); - //new TH2F("h2317","P2 vs X EXM cut st. 3",30,3.0,183.0,100,0.,1.5); - new TH2F("h2317","P2 vs X EXM cut st. 3",30,3.0,183.0,100,0.,6.); - new TH2F("h2318","P2 vs X EYM cut st. 3",30,3.0,183.0,100,0.,7.); - - // 2 - new TH1F("h2201","P2",30,3.0,183.0); - new TH2F("h2202","P2 vs X EPH2 st. 2",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2203","P2 vs X EAL2 st. 2",30,3.0,183.0,100,0.,0.005); - //new TH2F("h2204","P2 vs X EXM2 st. 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2204","P2 vs X EXM2 st. 2",30,3.0,183.0,100,0.,7.); - new TH2F("h2205","P2 vs X EYM2 st. 2",30,3.0,183.0,100,0.,5.); - - new TH2F("h2207","P vs X EPH2 st. 2",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2208","P vs X EAL2 st. 2",30,3.0,183.0,100,0.,0.005); - //new TH2F("h2209","P vs X EXM2 st. 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2209","P vs X EXM2 st. 2",30,3.0,183.0,100,0.,7.); - new TH2F("h2210","P vs X EYM2 st. 2",30,3.0,183.0,100,0.,5.); - - new TH2F("h2211","P vs X EPH cut st. 2",30,3.0,183.0,100,0.,0.05); - new TH2F("h2212","P vs X EAL cut st. 2",30,3.0,183.0,100,0.,0.2); - //new TH2F("h2213","P vs X EXM cut st. 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2213","P vs X EXM cut st. 2",30,3.0,183.0,100,0.,11.); - new TH2F("h2214","P vs X EYM cut st. 2",30,3.0,183.0,100,0.,10.); - - new TH2F("h2215","P2 vs X EPH cut st. 2",30,3.0,183.0,100,0.,0.05); - new TH2F("h2216","P2 vs X EAL cut st. 2",30,3.0,183.0,100,0.,0.2); - //new TH2F("h2217","P2 vs X EXM cut st. 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2217","P2 vs X EXM cut st. 2",30,3.0,183.0,100,0.,11.); - new TH2F("h2218","P2 vs X EYM cut st. 2",30,3.0,183.0,100,0.,10.); - - // 1 - new TH2F("h2102","P2 vs X EPH2 st. 2",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2103","P2 vs X EAL2 st. 2",30,3.0,183.0,100,0.,0.005); - //new TH2F("h2104","P2 vs X EXM2 st. 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2104","P2 vs X EXM2 st. 2",30,3.0,183.0,100,0.,7.); - new TH2F("h2105","P2 vs X EYM2 st. 2",30,3.0,183.0,100,0.,7.); - - new TH2F("h2107","P vs X EPH2 st. 2",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2108","P vs X EAL2 st. 2",30,3.0,183.0,100,0.,0.005); - //new TH2F("h2109","P vs X EXM2 st. 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2109","P vs X EXM2 st. 2",30,3.0,183.0,100,0.,7.); - new TH2F("h2110","P vs X EYM2 st. 2",30,3.0,183.0,100,0.,7.); - - new TH2F("h2111","P vs X EPH cut st. 2",30,3.0,183.0,100,0.,0.1); - new TH2F("h2112","P vs X EAL cut st. 2",30,3.0,183.0,100,0.,0.2); - //new TH2F("h2113","P vs X EXM cut st. 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2113","P vs X EXM cut st. 2",30,3.0,183.0,100,0.,11.); - new TH2F("h2114","P vs X EYM cut st. 2",30,3.0,183.0,100,0.,11.); - - new TH2F("h2115","P2 vs X EPH cut st. 2",30,3.0,183.0,100,0.,0.1); - new TH2F("h2116","P2 vs X EAL cut st. 2",30,3.0,183.0,100,0.,0.2); - //new TH2F("h2117","P2 vs X EXM cut st. 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2117","P2 vs X EXM cut st. 2",30,3.0,183.0,100,0.,11.); - new TH2F("h2118","P2 vs X EYM cut st. 2",30,3.0,183.0,100,0.,11.); - - // 2,3,4,5 - new TH1F("h2701","P2 fit 2",30,3.0,183.0); - new TH2F("h2702","P2 vs X EPH2 st. 1 fit 2",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2703","P2 vs X EAL2 st. 1 fit 2",30,3.0,183.0,100,0.,0.005); - // new TH2F("h2704","P2 vs X EXM2 st. 1 fit 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2704","P2 vs X EXM2 st. 1 fit 2",30,3.0,183.0,100,0.,2.); - new TH2F("h2705","P2 vs X EYM2 st. 1 fit 2",30,3.0,183.0,100,0.,3.); - - new TH2F("h2707","P vs X EPH2 st. 1 fit 2",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2708","P vs X EAL2 st. 1 fit 2",30,3.0,183.0,100,0.,0.005); - //new TH2F("h2709","P vs X EXM2 st. 1 fit 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2709","P vs X EXM2 st. 1 fit 2",30,3.0,183.0,100,0.,2.); - new TH2F("h2710","P vs X EYM2 st. 1 fit 2",30,3.0,183.0,100,0.,3.); - - new TH2F("h2711","P vs X EPH cut st. 1 fit 2",30,3.0,183.0,100,0.,0.07); - new TH2F("h2712","P vs X EAL cut st. 1 fit 2",30,3.0,183.0,100,0.,0.2); - //new TH2F("h2713","P vs X EXM cut st. 1 fit 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2713","P vs X EXM cut st. 1 fit 2",30,3.0,183.0,100,0.,6.); - new TH2F("h2714","P vs X EYM cut st. 1 fit 2",30,3.0,183.0,100,0.,7.); - - new TH2F("h2715","P2 vs X EPH cut st. 1 fit 2",30,3.0,183.0,100,0.,0.07); - new TH2F("h2716","P2 vs X EAL cut st. 1 fit 2",30,3.0,183.0,100,0.,0.2); - //new TH2F("h2717","P2 vs X EXM cut st. 1 fit 2",30,3.0,183.0,100,0.,1.5); - new TH2F("h2717","P2 vs X EXM cut st. 1 fit 2",30,3.0,183.0,100,0.,6.); - new TH2F("h2718","P2 vs X EYM cut st. 1 fit 2",30,3.0,183.0,100,0.,7.); - - // 1,3,4,5 - new TH1F("h2801","P2 fit 1",30,3.0,183.0); - new TH2F("h2802","P2 vs X EPH2 st. 2 fit 1",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2803","P2 vs X EAL2 st. 2 fit 1",30,3.0,183.0,100,0.,0.005); - //new TH2F("h2804","P2 vs X EXM2 st. 2 fit 1",30,3.0,183.0,100,0.,1.5); - new TH2F("h2804","P2 vs X EXM2 st. 2 fit 1",30,3.0,183.0,100,0.,2.); - new TH2F("h2805","P2 vs X EYM2 st. 2 fit 1",30,3.0,183.0,100,0.,3.); - - new TH2F("h2807","P vs X EPH2 st. 2 fit 1",30,3.0,183.0,100,0.,0.0006); - new TH2F("h2808","P vs X EAL2 st. 2 fit 1",30,3.0,183.0,100,0.,0.005); - //new TH2F("h2809","P vs X EXM2 st. 2 fit 1",30,3.0,183.0,100,0.,1.5); - new TH2F("h2809","P vs X EXM2 st. 2 fit 1",30,3.0,183.0,100,0.,2.); - new TH2F("h2810","P vs X EYM2 st. 2 fit 1",30,3.0,183.0,100,0.,3.); - - new TH2F("h2811","P vs X EPH cut st. 2 fit 1",30,3.0,183.0,100,0.,0.05); - new TH2F("h2812","P vs X EAL cut st. 2 fit 1",30,3.0,183.0,100,0.,0.2); - //new TH2F("h2813","P vs X EXM cut st. 2 fit 1",30,3.0,183.0,100,0.,1.5); - new TH2F("h2813","P vs X EXM cut st. 2 fit 1",30,3.0,183.0,100,0.,5.); - new TH2F("h2814","P vs X EYM cut st. 2 fit 1",30,3.0,183.0,100,0.,7.); - - new TH2F("h2815","P2 vs X EPH cut st. 2 fit 1",30,3.0,183.0,100,0.,0.05); - new TH2F("h2816","P2 vs X EAL cut st. 2 fit 1",30,3.0,183.0,100,0.,0.2); - //new TH2F("h2817","P2 vs X EXM cut st. 2 fit 1",30,3.0,183.0,100,0.,1.5); - new TH2F("h2817","P2 vs X EXM cut st. 2 fit 1",30,3.0,183.0,100,0.,5.); - new TH2F("h2818","P2 vs X EYM cut st. 2 fit 1",30,3.0,183.0,100,0.,7.); - // fin de test - - new TH1F("h500","Acceptance en H st. 4",500,0.,500.); - new TH1F("h600","Acceptance en H st. 5",500,0.,500.); - new TH1F("h700","X vertex track found",200,-10.,10.); - new TH1F("h701","Y vertex track found",200,-10.,10.); - new TH1F("h800","Rap. muon gen.",100,0.,5.); - new TH1F("h801","Rap. muon gen. recons.",100,0.,5.); - new TH1F("h802","Rap. muon gen. ghost ",100,0.,5.); - new TH1F("h900","Pt muon gen.",100,0.,20.); - new TH1F("h901","Pt muon gen. recons.",100,0.,20.); - new TH1F("h902","Pt muon gen. ghost",100,0.,20.); - new TH1F("h910","phi muon gen.",100,-10.,10.); - new TH1F("h911","phi muon gen. recons.",100,-10.,10.); - new TH1F("h912","phi muon gen. ghost",100,-10.,10.); - new TH2F("h1001","Y VS X hit st. 1",300,-300.,300.,300,-300.,300.); - new TH2F("h1002","Y VS X hit st. 2",300,-300.,300.,300,-300.,300.); - new TH2F("h1003","Y VS X hit st. 3",300,-300.,300.,300,-300.,300.); - new TH2F("h1004","Y VS X hit st. 4",300,-300.,300.,300,-300.,300.); - new TH2F("h1005","Y VS X hit st. 5",300,-300.,300.,300,-300.,300.); - // Histos variance dans 4 - new TH2F("h11","VAR X st. 4",30,3.0,183.0,100,0.,2.); - new TH2F("h12","VAR Y st. 4",30,3.0,183.0,100,0.,600.); - new TH2F("h13","VAR PHI st. 4",30,3.0,183.0,100,0.,0.0001); - new TH2F("h14","VAR ALM st. 4",30,3.0,183.0,100,0.,0.05); - new TH1F("h15","P",30,3.0,183.0); - new TH1F("h411","VAR X st. 4",100,-1.42,1.42); - new TH1F("h412","VAR Y st. 4",100,-25.,25.); - new TH1F("h413","VAR PHI st. 4",100,-0.01,0.01); - new TH1F("h414","VAR ALM st. 4",100,-0.23,0.23); - // histo2 - new TH2F("h211","histo2-VAR X st. 4",30,3.0,183.0,100,0.,2.); - new TH2F("h212","histo2-VAR Y st. 4",30,3.0,183.0,100,0.,600.); - new TH1F("h213","histo2-VAR X st. 4",100,-1.42,1.42); - new TH1F("h214","histo2-VAR Y st. 4",100,-25.,25.); - new TH1F("h215","histo2-P",30,3.0,183.0); - - // Histos variance dans 2 - new TH2F("h21","VAR X st. 2",30,3.0,183.0,100,0.,3.); - new TH2F("h22","VAR Y st. 2",30,3.0,183.0,100,0.,7.); - new TH2F("h23","VAR PHI st. 2",30,3.0,183.0,100,0.,0.006); - new TH2F("h24","VAR ALM st. 2",30,3.0,183.0,100,0.,0.005); - new TH1F("h25","P",30,3.0,183.0); - new TH1F("h421","VAR X st. 2",100,-1.72,1.72); - new TH1F("h422","VAR Y st. 2",100,-2.7,2.7); - new TH1F("h423","VAR PHI st. 2",100,-0.08,0.08); - new TH1F("h424","VAR ALM st. 2",100,-0.072,0.072); - // histo2 - new TH2F("h221","histo2-VAR X st. 2",30,3.0,183.0,100,0.,3.); - new TH2F("h222","histo2-VAR Y st. 2",30,3.0,183.0,100,0.,7.); - new TH1F("h223","histo2-VAR X st. 2",100,-1.72,1.72); - new TH1F("h224","histo2-VAR Y st. 2",100,-2.7,2.7); - new TH1F("h225","histo2-P",30,3.0,183.0); - - // Histos variance dans 1 - new TH2F("h31","VAR X st. 1",30,3.0,183.0,100,0.,2.); - new TH2F("h32","VAR Y st. 1",30,3.0,183.0,100,0.,0.5); - new TH2F("h33","VAR PHI st. 1",30,3.0,183.0,100,0.,0.006); - new TH2F("h34","VAR ALM st. 1",30,3.0,183.0,100,0.,0.005); - new TH1F("h35","P",30,3.0,183.0); - new TH1F("h431","VAR X st. 1",100,-1.42,1.42); - new TH1F("h432","VAR Y st. 1",100,-0.72,0.72); - new TH1F("h433","VAR PHI st. 1",100,-0.08,0.08); - new TH1F("h434","VAR ALM st. 1",100,-0.072,0.072); - // Histos variance dans 1 - new TH2F("h41","VAR X st. 1 fit 5,4,3,2,V",30,3.0,183.0,100,0.,4.); - new TH2F("h42","VAR Y st. 1 fit 5,4,3,2,V",30,3.0,183.0,100,0.,20.); - new TH2F("h43","VAR PHI st. 1 fit 5,4,3,2,V",30,3.0,183.0,100,0.,0.005); - new TH2F("h44","VAR ALM st. 1 fit 5,4,3,2,V",30,3.0,183.0,100,0.,0.005); - new TH1F("h45","P",30,3.0,183.0); - new TH1F("h441","VAR X st. 1 fit 5,4,3,2,V",100,-2.,2.); - new TH1F("h442","VAR Y st. 1 fit 5,4,3,2,V",100,-4.5,4.5); - new TH1F("h443","VAR PHI st. 1 fit 5,4,3,2,V",100,-0.072,0.072); - new TH1F("h444","VAR ALM st. 1 fit 5,4,3,2,V",100,-0.072,0.072); - // histo2 - new TH2F("h241","histo2-VAR X st. 1 fit 5,4,3,2,V",30,3.0,183.0,100,0.,4.); - new TH2F("h242","histo2-VAR Y st. 1 fit 5,4,3,2,V",30,3.0,183.0,100,0.,20.); - new TH1F("h243","histo2-VAR X st. 1 fit 5,4,3,2,V",100,-2.,2.); - new TH1F("h244","histo2-VAR Y st. 1 fit 5,4,3,2,V",100,-4.5,4.5); - new TH1F("h245","histo2-P",30,3.0,183.0); - - // Histos variance dans 2 - new TH2F("h51","VAR X st. 2 fit 5,4,3,1,V",30,3.0,183.0,100,0.,0.5); - new TH2F("h52","VAR Y st. 2 fit 5,4,3,1,V",30,3.0,183.0,100,0.,2.); - new TH2F("h53","VAR PHI st. 2 fit 5,4,3,1,V",30,3.0,183.0,100,0.,0.005); - new TH2F("h54","VAR ALM st. 2 fit 5,4,3,1,V",30,3.0,183.0,100,0.,0.01); - new TH1F("h55","P",30,3.0,183.0); - new TH1F("h451","VAR X st. 2 fit 5,4,3,1,V",100,-0.72,0.72); - new TH1F("h452","VAR Y st. 2 fit 5,4,3,1,V",100,-1.42,1.42); - new TH1F("h453","VAR PHI st. 2 fit 5,4,3,1,V",100,-0.072,0.072); - new TH1F("h454","VAR ALM st. 2 fit 5,4,3,1,V",100,-0.1,0.1); - new TH1F("h999","PTOT",30,3.0,183.0); - // histo2 - new TH2F("h251","histo2-VAR X st. 2 fit 5,4,3,1,V",30,3.0,183.0,100,0.,0.5); - new TH2F("h252","histo2-VAR Y st. 2 fit 5,4,3,1,V",30,3.0,183.0,100,0.,2.); - new TH1F("h253","histo2-VAR X st. 2 fit 5,4,3,1,V",100,-0.72,0.72); - new TH1F("h254","histo2-VAR Y st. 2 fit 5,4,3,1,V",100,-1.42,1.42); - new TH1F("h255","histo2-P",30,3.0,183.0); - // Histos variance dans 3 - new TH2F("h61","VAR X st. 3 fit 4,5,V",30,3.0,183.0,100,0.,5.); - new TH2F("h62","VAR Y st. 3 fit 4,5,V",30,3.0,183.0,100,0.,2.); - new TH2F("h63","VAR PHI st. 3 fit 4,5,V",30,3.0,183.0,100,0.,0.0006); - new TH2F("h64","VAR ALM st. 3 fit 4,5,V",30,3.0,183.0,100,0.,0.0006); - new TH1F("h65","P",30,3.0,183.0); - new TH1F("h461","VAR X st. 3 fit 4,5,V",100,-2.25,2.25); - new TH1F("h462","VAR Y st. 3 fit 4,5,V",100,-1.42,1.42); - new TH1F("h463","VAR PHI st. 3 fit 4,5,V",100,-0.024,0.024); - new TH1F("h464","VAR ALM st. 3 fit 4,5,V",100,-0.024,0.024); - // histo2 - new TH2F("h261","histo2-VAR X st. 3 fit 4,5,V",30,3.0,183.0,100,0.,5.); - new TH2F("h262","histo2-VAR Y st. 3 fit 4,5,V",30,3.0,183.0,100,0.,2.); - new TH1F("h263","histo2-VAR X st. 3 fit 4,5,V",100,-2.25,2.25); - new TH1F("h264","histo2-VAR Y st. 3 fit 4,5,V",100,-1.42,1.42); - new TH1F("h265","Phisto2-",30,3.0,183.0); - // Histos dx,dy distribution between chambers inside stations - new TH1F("h71","DX in st. ID-70",100,-5.,5.); - new TH1F("h81","DY in st. ID-80",100,-5.,5.); - new TH1F("h72","DX in st. ID-70",100,-5.,5.); - new TH1F("h82","DY in st. ID-80",100,-5.,5.); - new TH1F("h73","DX in st. ID-70",100,-5.,5.); - new TH1F("h83","DY in st. ID-80",100,-5.,5.); - new TH1F("h74","DX in st. ID-70",100,-5.,5.); - new TH1F("h84","DY in st. ID-80",100,-5.,5.); - new TH1F("h75","DX in st. ID-70",100,-5.,5.); - new TH1F("h85","DY in st. ID-80",100,-5.,5.); -} - -void chfnt(Int_t &ievr, Int_t &ntrackr, Int_t *istatr, Int_t *isignr, Float_t *pxr, Float_t *pyr, Float_t *pzr, Float_t *zvr, Float_t *chi2r, Float_t *pxv, Float_t *pyv, Float_t *pzv) -{ - // - // fill the ntuple - ntuple_st.ievr = ievr; - ntuple_st.ntrackr = ntrackr; - for (Int_t i=0; i<500; i++) { - ntuple_st.istatr[i] = istatr[i]; - ntuple_st.isignr[i] = isignr[i]; - ntuple_st.pxr[i] = pxr[i]; - ntuple_st.pyr[i] = pyr[i]; - ntuple_st.pzr[i] = pzr[i]; - ntuple_st.zvr[i] = zvr[i]; - ntuple_st.chi2r[i] = chi2r[i]; - ntuple_st.pxv[i] = pxv[i]; - ntuple_st.pyv[i] = pyv[i]; - ntuple_st.pzv[i] = pzv[i]; - } - ntuple_global->Fill(); -} - -void hist_closed() -{ - // - // write histos and ntuple to "reconst.root" file - hfile_global->Write(); -} - -void trackf_read_geant(Int_t *itypg, Double_t *xtrg, Double_t *ytrg, Double_t *ptotg, Int_t *idg, Int_t *izch, Double_t *pvert1g, Double_t *pvert2g, Double_t *pvert3g, Double_t *zvertg, Int_t &nhittot1, Double_t *cx, Double_t *cy, Double_t *cz, Int_t &ievr,Int_t &nev,Double_t *xgeant, Double_t *ygeant,Double_t *clsize1, Double_t *clsize2) -{ - // - // introduce aliroot variables in fortran common - // tracking study from geant hits - // - - AliMUON *MUON = (AliMUON*) gAlice->GetModule("MUON"); - - // TTree *TK = gAlice->TreeK(); - TTree *TH = gAlice->TreeH(); - Int_t ntracks = (Int_t)TH->GetEntries(); - cout<<"ntrack="<ResetHits(); - TH->GetEvent(track); - - if (MUON) { -// -// Loop over hits -// - for(AliMUONhit* mHit=(AliMUONhit*)MUON->FirstHit(-1); - mHit; - mHit=(AliMUONhit*)MUON->NextHit()) - { - if (maxidg<=20000) { - - if (mHit->fChamber > 10) continue; - TClonesArray *fPartArray = gAlice->Particles(); - TParticle *Part; - Int_t ftrack = mHit->fTrack; - Int_t id = ((TParticle*) fPartArray->UncheckedAt(ftrack))->GetPdgCode(); - - if (id==kMuonPlus||id==kMuonMinus) { - - // inversion de x et y car le champ est inverse dans le programme tracking - xtrg[maxidg] = 0; - ytrg[maxidg] = 0; - xgeant[maxidg] = mHit->fY; // x-pos of hit - ygeant[maxidg] = mHit->fX; // y-pos of hit - clsize1[maxidg] = 0; // cluster size on 1-st cathode - clsize2[maxidg] = 0; // cluster size on 2-nd cathode - cx[maxidg] = mHit->fCyHit; // Px/P of hit - cy[maxidg] = mHit->fCxHit; // Py/P of hit - cz[maxidg] = mHit->fCzHit; // Pz/P of hit - izch[maxidg] = mHit->fChamber; - /* - Int_t pdgtype = Int_t(mHit->fParticle); // particle number - itypg[maxidg] = gMC->IdFromPDG(pdgtype); - - */ - if (id==kMuonPlus) itypg[maxidg] = 5; - else itypg[maxidg] = 6; - - ptotg[maxidg] = mHit->fPTot; // P of hit - - Part = (TParticle*) fPartArray->UncheckedAt(ftrack); - Float_t thet = Part->Theta(); - thet = thet*180./3.1416; - - Int_t iparent = Part->GetFirstMother(); - if (iparent >= 0) { - Int_t ip; - while(1) { - ip=((TParticle*) fPartArray->UncheckedAt(iparent))->GetFirstMother(); - if (ip < 0) { - break; - } else { - iparent = ip; - } - } - } - //printf("iparent - %d\n",iparent); - Int_t id1 = ftrack; // numero de la particule generee au vertex - Int_t idum = track+1; - Int_t id2 = ((TParticle*) fPartArray->UncheckedAt(iparent))->GetPdgCode(); - - if (id2==443) id2=114; - else id2=116; - - if (id2==116) { - nres++; - } - //printf("id2 %d\n",id2); - idg[maxidg] = 30000*id1+10000*idum+id2; - - pvert1g[maxidg] = Part->Py(); // Px vertex - pvert2g[maxidg] = Part->Px(); // Py vertex - pvert3g[maxidg] = Part->Pz(); // Pz vertex - zvertg[maxidg] = Part->Vz(); // z vertex - maxidg ++; - - } - } - } // hit loop - } // if MUON - } // track loop first file - - if (TrH1 && fHits2 ) { // if background file - ntracks =(Int_t)TrH1->GetEntries(); - printf("Trackf_read - 2-nd file - ntracks %d\n",ntracks); - - // Loop over tracks - for (Int_t track=0; trackClear(); - TrH1->GetEvent(track); - - // Loop over hits - for (int i=0;iGetEntriesFast();i++) - { - AliMUONhit *mHit=(AliMUONhit*) (*fHits2)[i]; - - if (mHit->fChamber > 10) continue; - - if (maxidg<=20000) { - - // inversion de x et y car le champ est inverse dans le programme tracking !!!! - xtrg[maxidg] = 0; // only for reconstructed point - ytrg[maxidg] = 0; // only for reconstructed point - xgeant[maxidg] = mHit->fY; // x-pos of hit - ygeant[maxidg] = mHit->fX; // y-pos of hit - clsize1[maxidg] = 0; // cluster size on 1-st cathode - clsize2[maxidg] = 0; // cluster size on 2-nd cathode - cx[maxidg] = mHit->fCyHit; // Px/P of hit - cy[maxidg] = mHit->fCxHit; // Py/P of hit - cz[maxidg] = mHit->fCzHit; // Pz/P of hit - izch[maxidg] = mHit->fChamber; // chamber number - ptotg[maxidg] = mHit->fPTot; // P of hit - - Int_t ftrack = mHit->fTrack; - Int_t id1 = ftrack; // track number - Int_t idum = track+1; - - TClonesArray *fPartArray = fParticles2; - TParticle *Part; - Part = (TParticle*) fPartArray->UncheckedAt(ftrack); - Int_t id = ((TParticle*) fPartArray->UncheckedAt(ftrack))->GetPdgCode(); - if (id==kMuonPlus||id==kMuonMinus) { - if (id==kMuonPlus) itypg[maxidg] = 5; - else itypg[maxidg] = 6; - } else itypg[maxidg]=0; - - Int_t id2=0; // set parent to 0 for background !! - idg[maxidg] = 30000*id1+10000*idum+id2; - - pvert1g[maxidg] = Part->Py(); // Px vertex - pvert2g[maxidg] = Part->Px(); // Py vertex - pvert3g[maxidg] = Part->Pz(); // Pz vertex - zvertg[maxidg] = Part->Vz(); // z vertex - - maxidg ++; - - } // check limits (maxidg) - } // hit loop - } // track loop - } // if TrH1 - - ievr = nev; - nhittot1 = maxidg ; - cout<<"nhittot1="<=19) nbres++; - printf("nres, nbres %d %d \n",nres,nbres); - - hfile_global->cd(); - -} - - - -void trackf_read_spoint(Int_t *itypg, Double_t *xtrg, Double_t *ytrg, Double_t *ptotg, Int_t *idg, Int_t *izch, Double_t *pvert1g, Double_t *pvert2g, Double_t *pvert3g, Double_t *zvertg, Int_t &nhittot1, Double_t *cx, Double_t *cy, Double_t *cz, Int_t &ievr,Int_t &nev,Double_t *xgeant, Double_t *ygeant,Double_t *clsize1, Double_t *clsize2) - -{ - // - // introduce aliroot variables in fortran common - // tracking study from reconstructed points - // - AliMUON *MUON = (AliMUON*) gAlice->GetModule("MUON"); - - cout<<"numero de l'evenement "<GetTreeC(nev); - TTree *TC=MUON->TreeC(); - TC->GetEntries(); - - Int_t maxidg = 0; - Int_t nres=0; - Int_t nncor=0; - static Int_t nuncor=0; - static Int_t nbadcor=0; - AliMUONRawCluster * mRaw; - AliMUONRawCluster * mRaw1; - TTree *TH = gAlice->TreeH(); - - Int_t ihit; - Int_t mult1, mult2; - if (MUON) { - for (Int_t ich=0;ich<10;ich++) { - TClonesArray *MUONcorrel = MUON->CathCorrelAddress(ich); - MUON->ResetCorrelation(); - TC->GetEvent(); - Int_t ncor = (Int_t)MUONcorrel->GetEntries(); - if (ncor>=2) nncor++; - if (!ncor) continue; - - // Loop over correlated clusters - for (Int_t icor=0;icorUncheckedAt(icor); - - Int_t flag=0; // = 1 if no information in the second cathode - Int_t index = mCor->fCorrelIndex[0]; // for the second cathode - if (index >= 0) { - Int_t index1 = mCor->fCorrelIndex[3]; // for the 1-st cathode - mRaw1 = MUON->RawCluster(ich,1,index1); - mult1=mRaw1->fMultiplicity; - mRaw = MUON->RawCluster(ich,2,index); - mult2=mRaw->fMultiplicity; - } else { - index = mCor->fCorrelIndex[3]; - mRaw = MUON->RawCluster(ich,1,index); - mult1=mRaw->fMultiplicity; - mult2=0; - flag=1; - nuncor++; - } - if (!mRaw) continue; - - Int_t ftrack1 = mRaw->fTracks[1]; // qui doit etre le meme pour - // la cathode 1 et 2 - ihit= mRaw->fTracks[0]; - //printf("icor, ftrack1 ihit %d %d %d\n",icor,ftrack1,ihit); - - if (mRaw->fClusterType == 0 ) { - - if (maxidg<=20000) { - if (flag == 0) { - xtrg[maxidg] = (Double_t) mCor->fY[3]; - ytrg[maxidg] = (Double_t) mCor->fX[0]; - Int_t index1 = mCor->fCorrelIndex[3]; - mRaw1 = MUON->RawCluster(ich,1,index1); - if (mRaw1->fClusterType==1 || mRaw1->fClusterType==2) { - Float_t xclust=mCor->fX[3]; - Float_t yclust=mCor->fY[3]; - AliMUONchamber *iChamber=&(MUON->Chamber(ich)); - AliMUONsegmentation *seg = iChamber->GetSegmentationModel(1); - Int_t ix,iy; - seg->GetPadIxy(xclust,yclust,ix,iy); - Int_t isec=seg->Sector(ix,iy); - printf("nev, CORRELATION with pure background in chamber sector %d %d %d !!!!!!!!!!!!!!!!!!!!!\n",nev,ich+1,isec); - nbadcor++; - - } // end if cluster type on cathode 1 - }else { - xtrg[maxidg] = (Double_t) mCor->fY[3]; - ytrg[maxidg] = (Double_t) mCor->fX[3]; - } // if iflag - izch[maxidg] = ich+1; - xgeant[maxidg] = 0; - ygeant[maxidg] = 0; - clsize1[maxidg] = mult1; - clsize2[maxidg] = mult2; - - cx[maxidg] = 0; // Px/P of hit - cy[maxidg] = 0; // Py/P of hit - cz[maxidg] = 0; // Pz/P of hit - itypg[maxidg] = 0; // particle number - ptotg[maxidg] = 0; // P of hit - idg[maxidg] = 0; - pvert1g[maxidg] = 0; // Px vertex - pvert2g[maxidg] = 0; // Py vertex - pvert3g[maxidg] = 0; // Pz vertex - zvertg[maxidg] = 0; // z vertex - maxidg++; - - }// fin maxidg - - } else if (mRaw->fClusterType ==1 && ftrack1 < 0) // background + resonance - { - nres++; - // get indexmap and loop over digits to find the signal - Int_t nent=(Int_t)gAlice->TreeD()->GetEntries(); - gAlice->ResetDigits(); - if (flag==0) { - //gAlice->TreeD()->GetEvent(2); // cathode 2 - gAlice->TreeD()->GetEvent(nent-1); // cathode 2 - } else { - //gAlice->TreeD()->GetEvent(1); // cathode 1 - gAlice->TreeD()->GetEvent(nent-2); // cathode 1 - } - - TClonesArray *MUONdigits = MUON->DigitsAddress(ich); - Int_t mul=mRaw->fMultiplicity; - Int_t trsign; - for (int i=0;ifIndexMap[i]; - AliMUONdigit *dig= (AliMUONdigit*)MUONdigits->UncheckedAt(idx); - trsign=dig->fTracks[0]; - ihit=dig->fHit-1; - if (trsign > 0 && ihit >= 0) break; - - } // loop over indexmap - - //printf("trsign, ihit %d %d\n",trsign,ihit); - //printf("signal+background : trsign %d\n",trsign); - - if (trsign < 0 || ihit < 0) { // no signal muon was found - - if (maxidg<=20000) { - if (flag == 0) { - xtrg[maxidg] = (Double_t) mCor->fY[3]; - ytrg[maxidg] = (Double_t) mCor->fX[0]; - }else { - xtrg[maxidg] = (Double_t) mCor->fY[3]; - ytrg[maxidg] = (Double_t) mCor->fX[3]; - } - - izch[maxidg] = ich+1; - - // initialisation of informations which - // can't be reached for background - - xgeant[maxidg] = 0; // only for resonances - ygeant[maxidg] = 0; // only for resonances - clsize1[maxidg] = mult1; - clsize2[maxidg] = mult2; - - cx[maxidg] = 0; // Px/P of hit - cy[maxidg] = 0; // Py/P of hit - cz[maxidg] = 0; // Pz/P of hit - itypg[maxidg] = 0; // particle number - ptotg[maxidg] = 0; // P of hit - idg[maxidg] = 0; - pvert1g[maxidg] = 0; // Px vertex - pvert2g[maxidg] = 0; // Py vertex - pvert3g[maxidg] = 0; // Pz vertex - zvertg[maxidg] = 0; - maxidg++; - - }// fin maxidg - } else { // signal muon - retrieve info - //printf("inside trsign, ihit %d %d\n",trsign,ihit); - if (maxidg<=20000) { - if (flag == 0) { - xtrg[maxidg] = (Double_t) mCor->fY[3]; - ytrg[maxidg] = (Double_t) mCor->fX[0]; - }else { - xtrg[maxidg] = (Double_t) mCor->fY[3]; - ytrg[maxidg] = (Double_t) mCor->fX[3]; - } - izch[maxidg] = ich+1; - clsize1[maxidg] = mult1; - clsize2[maxidg] = mult2; - - // initialise and set to the correct values - // if signal muons - - xgeant[maxidg] = 0; // only for resonances - ygeant[maxidg] = 0; // only for resonances - - cx[maxidg] = 0; // Px/P of hit - cy[maxidg] = 0; // Py/P of hit - cz[maxidg] = 0; // Pz/P of hit - itypg[maxidg] = 0; // particle number - ptotg[maxidg] = 0; // P of hit - idg[maxidg] = 0; - pvert1g[maxidg] = 0; // Px vertex - pvert2g[maxidg] = 0; // Py vertex - pvert3g[maxidg] = 0; // Pz vertex - zvertg[maxidg] = 0; - // try to retrieve info about signal muons - gAlice->ResetHits(); - TH->GetEvent(trsign); - - TClonesArray *MUONhits = MUON->Hits(); - AliMUONhit *mHit= (AliMUONhit*)MUONhits-> - UncheckedAt(ihit); - TClonesArray *fPartArray = gAlice->Particles(); - TParticle *Part; - Int_t nch=mHit->fChamber-1; - //printf("sig+bgr ich, nch %d %d \n",ich,nch); - if (nch==ich) { - Int_t ftrack = mHit->fTrack; - Int_t id = ((TParticle*) fPartArray-> - UncheckedAt(ftrack))->GetPdgCode(); - if (id==kMuonPlus||id==kMuonMinus) { - xgeant[maxidg] = (Double_t) mHit->fY; - ygeant[maxidg] = (Double_t) mHit->fX; - cx[maxidg] = (Double_t) mHit->fCyHit; - cy[maxidg] = (Double_t) mHit->fCxHit; - cz[maxidg] = (Double_t) mHit->fCzHit; - - if (id==kMuonPlus) { - itypg[maxidg] = 5; - } else if (id==kMuonMinus) { - itypg[maxidg] = 6; - } else itypg[maxidg] = 0; - - ptotg[maxidg] = (Double_t) mHit->fPTot; - Part = (TParticle*) fPartArray-> - UncheckedAt(ftrack); - Int_t iparent = Part->GetFirstMother(); - Int_t id2; - id2 = ((TParticle*) fPartArray-> - UncheckedAt(ftrack))->GetPdgCode(); - - if (iparent >= 0) { - Int_t ip; - while(1) { - ip=((TParticle*) fPartArray-> - UncheckedAt(iparent))->GetFirstMother(); - if (ip < 0) { - id2 = ((TParticle*) fPartArray-> - UncheckedAt(iparent))->GetPdgCode(); - break; - } else { - iparent = ip; - id2 = ((TParticle*) fPartArray-> - UncheckedAt(iparent))->GetPdgCode(); - } // ip<0 - } // while - }// iparent - Int_t id1 = ftrack; - Int_t idum = trsign+1; - - if (id2==443 || id2==553) { - nres++; - if (id2==443) id2=114; - else id2=116; - } - - idg[maxidg] = 30000*id1+10000*idum+id2; - pvert1g[maxidg] = (Double_t) Part->Py(); - pvert2g[maxidg] = (Double_t) Part->Px(); - pvert3g[maxidg] = (Double_t) Part->Pz(); - zvertg[maxidg] = (Double_t) Part->Vz(); - } //if muon - } //if nch - maxidg++; - } // check limits - } // sign+bgr, highest bgr - } - //pure resonance or mixed cluster with the highest - //contribution coming from resonance - if (mRaw->fClusterType >= 1 && ftrack1>=0) - { - if (maxidg<=20000) { - if (flag == 0) { - xtrg[maxidg] = (Double_t) mCor->fY[3]; - ytrg[maxidg] = (Double_t) mCor->fX[0]; - }else { - xtrg[maxidg] = (Double_t) mCor->fY[3]; - ytrg[maxidg] = (Double_t) mCor->fX[3]; - } - clsize1[maxidg] = mult1; - clsize2[maxidg] = mult2; - izch[maxidg] = ich+1; - - Int_t nent=(Int_t)gAlice->TreeD()->GetEntries(); - gAlice->ResetDigits(); - if (flag==0) { - //gAlice->TreeD()->GetEvent(2); // cathode 2 - gAlice->TreeD()->GetEvent(nent-1); // cathode 2 - } else { - //gAlice->TreeD()->GetEvent(1); // cathode 1 - gAlice->TreeD()->GetEvent(nent-2); // cathode 1 - } - - TClonesArray *MUONdigits = MUON->DigitsAddress(ich); - Int_t mul=mRaw->fMultiplicity; - for (int i=0;ifIndexMap[i]; - AliMUONdigit *dig= (AliMUONdigit*)MUONdigits->UncheckedAt(idx); - ihit=dig->fHit-1; - if (ihit >= 0) break; - - } // loop over indexmap - //printf("fClusterType, ihit %d %d \n",mRaw->fClusterType,ihit); - if (ihit < 0) { - xgeant[maxidg] = 0; // only for resonances - ygeant[maxidg] = 0; // only for resonances - - cx[maxidg] = 0; // Px/P of hit - cy[maxidg] = 0; // Py/P of hit - cz[maxidg] = 0; // Pz/P of hit - itypg[maxidg] = 0; // particle number - ptotg[maxidg] = 0; // P of hit - idg[maxidg] = 0; - pvert1g[maxidg] = 0; // Px vertex - pvert2g[maxidg] = 0; // Py vertex - pvert3g[maxidg] = 0; // Pz vertex - zvertg[maxidg] = 0; - } else { - gAlice->ResetHits(); - TH->GetEvent(ftrack1); - TClonesArray *MUONhits = MUON->Hits(); - AliMUONhit *mHit= (AliMUONhit*)MUONhits-> - UncheckedAt(ihit); - TClonesArray *fPartArray = gAlice->Particles(); - TParticle *Part; - Int_t nch=mHit->fChamber-1; - //printf("signal ich, nch %d %d \n",ich,nch); - if (nch==ich) { - Int_t ftrack = mHit->fTrack; - Int_t id = ((TParticle*) fPartArray-> - UncheckedAt(ftrack))->GetPdgCode(); - //printf("id %d \n",id); - if (id==kMuonPlus||id==kMuonMinus) { - xgeant[maxidg] = (Double_t) mHit->fY; - ygeant[maxidg] = (Double_t) mHit->fX; - cx[maxidg] = (Double_t) mHit->fCyHit; - cy[maxidg] = (Double_t) mHit->fCxHit; - cz[maxidg] = (Double_t) mHit->fCzHit; - - if (id==kMuonPlus) { - itypg[maxidg] = 5; - } else if (id==kMuonMinus) { - itypg[maxidg] = 6; - } else itypg[maxidg] = 0; - - ptotg[maxidg] = (Double_t) mHit->fPTot; - Part = (TParticle*) fPartArray-> - UncheckedAt(ftrack); - Int_t iparent = Part->GetFirstMother(); - Int_t id2; - id2 = ((TParticle*) fPartArray-> - UncheckedAt(ftrack))->GetPdgCode(); - - if (iparent >= 0) { - Int_t ip; - while(1) { - ip=((TParticle*) fPartArray-> - UncheckedAt(iparent))->GetFirstMother(); - if (ip < 0) { - id2 = ((TParticle*) fPartArray-> - UncheckedAt(iparent))->GetPdgCode(); - break; - } else { - iparent = ip; - id2 = ((TParticle*) fPartArray-> - UncheckedAt(iparent))->GetPdgCode(); - } // ip<0 - } // while - }// iparent - Int_t id1 = ftrack; - Int_t idum = ftrack1+1; - - if (id2==443 || id2==553) { - nres++; - if (id2==443) id2=114; - else id2=116; - } - // printf("id2 %d\n",id2); - idg[maxidg] = 30000*id1+10000*idum+id2; - pvert1g[maxidg] = (Double_t) Part->Py(); - pvert2g[maxidg] = (Double_t) Part->Px(); - pvert3g[maxidg] = (Double_t) Part->Pz(); - zvertg[maxidg] = (Double_t) Part->Vz(); - } //if muon - } //if nch - } // ihit - maxidg++; - } // check limits - } // if cluster type - } // icor loop - } // ich loop - }// if MUON - - - ievr = nev; - cout<<"evenement "<=20) nbcor++; - printf("nbcor - %d\n",nbcor); - printf("nuncor - %d\n",nuncor); - printf("nbadcor - %d\n",nbadcor); - - TC->Reset(); - - hfile_global->cd(); - -} - -void trackf_fit(Int_t &ivertex, Double_t *pest, Double_t *pstep, Double_t &pxzinv, Double_t &tphi, Double_t &talam, Double_t &xvert, Double_t &yvert) -{ - // - // Fit a track candidate with the following input parameters: - // INPUT : IVERTEX : vertex flag, if IVERTEX=1 (XVERT,YVERT) are free paramaters - // if IVERTEX=1 (XVERT,YVERT)=(0.,0.) - // PEST(5) : starting value of parameters (minuit) - // PSTEP(5) : step size for parameters (minuit) - // OUTPUT : PXZINV,TPHI,TALAM,XVERT,YVERT : fitted value of the parameters - - static Double_t arglist[10]; - static Double_t c[5] = {0.4, 0.45, 0.45, 90., 90.}; - static Double_t b1, b2, epxz, efi, exs, exvert, eyvert; - TString chname; - Int_t ierflg = 0; - - TMinuit *gMinuit = new TMinuit(5); - gMinuit->mninit(5,10,7); - gMinuit->SetFCN(fcnfwrap); // constant m.f. - - arglist[0] = -1; - - gMinuit->mnexcm("SET PRINT", arglist, 1, ierflg); - // gMinuit->mnseti('track fitting'); - - gMinuit->mnparm(0, "invmom", pest[0], pstep[0], -c[0], c[0], ierflg); - gMinuit->mnparm(1, "azimuth", pest[1], pstep[1], -c[1], c[1], ierflg); - gMinuit->mnparm(2, "deep", pest[2], pstep[2], -c[2], c[2], ierflg); - if (ivertex==1) { - gMinuit->mnparm(3, "x ", pest[3], pstep[3], -c[3], c[3], ierflg); - gMinuit->mnparm(4, "y ", pest[4], pstep[4], -c[4], c[4], ierflg); - } - - gMinuit->mnexcm("SET NOGR", arglist, 0, ierflg); - gMinuit->mnexcm("MINIMIZE", arglist, 0, ierflg); - gMinuit->mnexcm("EXIT" , arglist, 0, ierflg); - - gMinuit->mnpout(0, chname, pxzinv, epxz, b1, b2, ierflg); - gMinuit->mnpout(1, chname, tphi, efi, b1, b2, ierflg); - gMinuit->mnpout(2, chname, talam, exs, b1, b2, ierflg); - if (ivertex==1) { - gMinuit->mnpout(3, chname, xvert, exvert, b1, b2, ierflg); - gMinuit->mnpout(4, chname, yvert, eyvert, b1, b2, ierflg); - } - - delete gMinuit; - -} - -void fcnf(Int_t &npar, Double_t *grad, Double_t &fval, Double_t *pest, Int_t iflag) -{ - // - // function called by trackf_fit - Int_t futil = 0; - fcn(npar,grad,fval,pest,iflag,futil); -} - -void prec_fit(Double_t &pxzinv, Double_t &fis, Double_t &alams, Double_t &xvert, Double_t &yvert, Double_t &pxzinvf, Double_t &fif, Double_t &alf, Double_t &xvertf, Double_t &yvertf, Double_t &epxzinv, Double_t &efi, Double_t &exs, Double_t &exvert, Double_t &eyvert) -{ - // - // minuit fits for tracking finding - - static Double_t arglist[10]; - static Double_t c1[5] = {0.001, 0.001, 0.001, 1., 1.}; - static Double_t c2[5] = {0.5, 0.5, 0.5, 120., 120.}; - static Double_t emat[9]; - static Double_t b1, b2; - Double_t fmin, fedm, errdef; - Int_t npari, nparx, istat; - - TString chname; - Int_t ierflg = 0; - - TMinuit *gMinuit = new TMinuit(5); - gMinuit->mninit(5,10,7); - gMinuit->SetFCN(fcnfitfwrap); - - arglist[0] = -1.; - gMinuit->mnexcm("SET PRINT", arglist, 1, ierflg); - - // gMinuit->mnseti('track fitting'); - - gMinuit->mnparm(0,"invmom", pxzinv, c1[0], -c2[0], c2[0], ierflg); // 0.003, 0.5 - gMinuit->mnparm(1,"azimuth ", fis, c1[1], -c2[1], c2[1], ierflg); - gMinuit->mnparm(2,"deep ", alams, c1[2], -c2[2], c2[2], ierflg); - gMinuit->mnparm(3,"xvert", xvert, c1[3], -c2[3], c2[3], ierflg); - gMinuit->mnparm(4,"yvert", yvert, c1[4], -c2[4], c2[4], ierflg); - - gMinuit->mnexcm("SET NOGR", arglist, 0, ierflg); - arglist[0] = 2.; - gMinuit->mnexcm("MINIMIZE", arglist, 0, ierflg); - gMinuit->mnexcm("EXIT", arglist, 0, ierflg); - - gMinuit->mnpout(0, chname, pxzinvf, epxzinv, b1, b2, ierflg); - gMinuit->mnpout(1, chname, fif, efi, b1, b2, ierflg); - gMinuit->mnpout(2, chname, alf, exs, b1, b2, ierflg); - gMinuit->mnpout(3, chname, xvertf, exvert, b1, b2, ierflg); - gMinuit->mnpout(4, chname, yvertf, eyvert, b1, b2, ierflg); - - gMinuit->mnemat(emat, 3); - gMinuit->mnstat(fmin, fedm, errdef, npari, nparx, istat); - - delete gMinuit; -} - -void fcnfitf(Int_t &npar, Double_t *grad, Double_t &fval, Double_t *xval, Int_t iflag) -{ - // - // function called by prec_fit - Int_t futil = 0; - fcnfit(npar,grad,fval,xval,iflag,futil); -} - -///////////////////// fin modifs perso ////////////////////// - -ClassImp(AliMUONcluster) - -//___________________________________________ -AliMUONcluster::AliMUONcluster(Int_t *clhits) -{ - fHitNumber=clhits[0]; - fCathode=clhits[1]; - fQ=clhits[2]; - fPadX=clhits[3]; - fPadY=clhits[4]; - fQpad=clhits[5]; - fRSec=clhits[6]; -} -ClassImp(AliMUONdigit) -//_____________________________________________________________________________ -AliMUONdigit::AliMUONdigit(Int_t *digits) -{ - // - // Creates a MUON digit object to be updated - // - fPadX = digits[0]; - fPadY = digits[1]; - fSignal = digits[2]; - fPhysics = digits[3]; - fHit = digits[4]; - -} -//_____________________________________________________________________________ -AliMUONdigit::AliMUONdigit(Int_t *tracks, Int_t *charges, Int_t *digits) -{ - // - // Creates a MUON digit object - // - fPadX = digits[0]; - fPadY = digits[1]; - fSignal = digits[2]; - fPhysics = digits[3]; - fHit = digits[4]; - for(Int_t i=0; i<10; i++) { - fTcharges[i] = charges[i]; - fTracks[i] = tracks[i]; - } -} - -AliMUONdigit::~AliMUONdigit() -{ - -} - -ClassImp(AliMUONlist) - -//____________________________________________________________________________ - AliMUONlist::AliMUONlist(Int_t ich, Int_t *digits): - AliMUONdigit(digits) +AliMUON& AliMUON::operator = (const AliMUON& rhs) { - // - // Creates a MUON digit list object - // - - fChamber = ich; - fTrackList = new TObjArray; - -} - -ClassImp(AliMUONhit) - -//___________________________________________ - AliMUONhit::AliMUONhit(Int_t shunt, Int_t track, Int_t *vol, Float_t *hits): - AliHit(shunt, track) -{ - fChamber=vol[0]; - fParticle=hits[0]; - fX=hits[1]; - fY=hits[2]; - fZ=hits[3]; - fTheta=hits[4]; - fPhi=hits[5]; - fTlength=hits[6]; - fEloss=hits[7]; - fPHfirst=(Int_t) hits[8]; - fPHlast=(Int_t) hits[9]; - - // modifs perso - fPTot=hits[10]; - fCxHit=hits[11]; - fCyHit=hits[12]; - fCzHit=hits[13]; -} -ClassImp(AliMUONcorrelation) -//___________________________________________ -//_____________________________________________________________________________ -AliMUONcorrelation::AliMUONcorrelation(Int_t *idx, Float_t *x, Float_t *y) -{ - // - // Creates a MUON correlation object - // - for(Int_t i=0; i<4; i++) { - fCorrelIndex[i] = idx[i]; - fX[i] = x[i]; - fY[i] = y[i]; - } +// copy operator +// dummy version + return *this; } -ClassImp(AliMUONRawCluster) -Int_t AliMUONRawCluster::Compare(TObject *obj) -{ - /* - AliMUONRawCluster *raw=(AliMUONRawCluster *)obj; - Float_t r=GetRadius(); - Float_t ro=raw->GetRadius(); - if (r>ro) return 1; - else if (rfY; - if (y>yo) return 1; - else if (y1) { - half=(high+low)/2; - if(y>coord[half]) low=half; - else high=half; - } - return low; -} - -void AliMUONRawCluster::SortMin(Int_t *idx,Float_t *xdarray,Float_t *xarray,Float_t *yarray,Float_t *qarray, Int_t ntr) -{ - // - // Get the 3 closest points(cog) one can find on the second cathode - // starting from a given cog on first cathode - // - - // - // Loop over deltax, only 3 times - // - - Float_t xmin; - Int_t jmin; - Int_t id[3] = {-2,-2,-2}; - Float_t jx[3] = {0.,0.,0.}; - Float_t jy[3] = {0.,0.,0.}; - Float_t jq[3] = {0.,0.,0.}; - Int_t jid[3] = {-2,-2,-2}; - Int_t i,j,imax; - - if (ntr<3) imax=ntr; - else imax=3; - for(i=0;i