New version from B. Batyunya to get the Dubna model work with the present HEAD
authorbarbera <barbera@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 14 Jun 2001 14:33:53 +0000 (14:33 +0000)
committerbarbera <barbera@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 14 Jun 2001 14:33:53 +0000 (14:33 +0000)
ITS/AliITSClusterFinderSPDdubna.cxx
ITS/AliITSRawCluster.cxx
ITS/AliITSsimulationSPDdubna.cxx
ITS/SPDclusterTestDubna.C

index 9d3d99c..0153f3d 100644 (file)
@@ -46,7 +46,7 @@ AliITSClusterFinderSPDdubna::AliITSClusterFinderSPDdubna
     fSegmentation=seg;
     fDigits=digits;
     fClusters=recp;
-    fNclusters= 0;
+    fNclusters=0;
     fMap=new AliITSMapA1(fSegmentation,fDigits);
     SetDx();
     SetDz();
@@ -74,8 +74,8 @@ AliITSClusterFinderSPDdubna::~AliITSClusterFinderSPDdubna()
   // destructor
   if (fMap) delete fMap;
 
-
 }
+
 //_____________________________________________________________________________
 
 void AliITSClusterFinderSPDdubna::Find1DClusters()
@@ -199,7 +199,6 @@ void  AliITSClusterFinderSPDdubna::GroupClusters()
   // get number of clusters for this module
   Int_t nofClusters = fClusters->GetEntriesFast();
   nofClusters -= fNclusters;
-
   AliITSRawClusterSPD *clusterI;
   AliITSRawClusterSPD *clusterJ;
   
@@ -239,6 +238,7 @@ void  AliITSClusterFinderSPDdubna::GroupClusters()
   } // I clusters
   fClusters->Compress();
 
+
   /*
     Int_t totalNofClusters = fClusters->GetEntriesFast();
     cout << " Nomber of clusters at the group end ="<< totalNofClusters<<endl;
@@ -361,18 +361,30 @@ void AliITSClusterFinderSPDdubna::GetRecPoints()
 
   Float_t spdLength = fSegmentation->Dz();
   Float_t spdWidth = fSegmentation->Dx();
-
+  Int_t dummy = 0;
+  Float_t xpitch = fSegmentation->Dpx(dummy);
   Int_t i;
   Int_t track0, track1, track2;
 
   for(i=0; i<nofClusters; i++) { 
 
     AliITSRawClusterSPD *clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
+    Int_t clustersizex = clusterI->NclX();
+    Int_t clustersizez = clusterI->NclZ();
+    Int_t xstart = clusterI->XStartf();
+    Float_t clusterx = 0;
+    Float_t clusterz = clusterI->Z();
+
+    for(Int_t ii=0; ii<clustersizex; ii++) {
+      clusterx += (xstart+0.5+ii)*xpitch;
+    }
+      clusterx /= clustersizex;
+      clusterz /= clustersizez;   
     clusterI->GetTracks(track0, track1, track2); 
     AliITSRecPoint rnew;
 
-    rnew.SetX((clusterI->X() - spdWidth/2)*kconv);
-    rnew.SetZ((clusterI->Z() - spdLength/2)*kconv);
+    rnew.SetX((clusterx - spdWidth/2)*kconv);
+    rnew.SetZ((clusterz - spdLength/2)*kconv);
     rnew.SetQ(1.);
     rnew.SetdEdX(0.);
     rnew.SetSigmaX2(kRMSx*kRMSx);
@@ -391,6 +403,7 @@ void AliITSClusterFinderSPDdubna::GetRecPoints()
 void AliITSClusterFinderSPDdubna::FindRawClusters(Int_t mod)
 {
   // find raw clusters
+  cout<<"SPDdubna: module ="<<mod<<endl;
   Find1DClusters();
   GroupClusters();
   TracksInCluster();
index a9e7f79..3b9e14a 100644 (file)
@@ -153,7 +153,7 @@ void AliITSRawClusterSPD::Add(AliITSRawClusterSPD* clJ) {
   
   if(this->fZStop < clJ->ZStop()) this->fZStop = clJ->ZStop();  
   
-  this->fZ = (this->fZ + clJ->Z())/2.;
+  this->fZ = this->fZ + clJ->Z();
   this->fX = (this->fX + clJ->X())/2.;
   this->fQ = this->fQ + clJ->Q();
   
index dd25d84..3326d96 100644 (file)
@@ -115,9 +115,10 @@ void AliITSsimulationSPDdubna::DigitiseModule(AliITSmodule *mod, Int_t module, I
 
     Float_t spdLength = fSegmentation->Dz();
     Float_t spdWidth = fSegmentation->Dx();
-
+    Float_t spdThickness = fSegmentation->Dy();
     Float_t difCoef, dum;       
     fResponse->DiffCoeff(difCoef,dum); 
+    if(spdThickness > 290) difCoef = 0.00613;  
 
     Float_t zPix0 = 1e+6;
     Float_t xPix0 = 1e+6;
@@ -130,7 +131,7 @@ void AliITSsimulationSPDdubna::DigitiseModule(AliITSmodule *mod, Int_t module, I
     Int_t nhits = fHits->GetEntriesFast();
     if (!nhits) return;
 
-    //cout<<"len,wid,dy,nx,nz,pitchx,pitchz ="<<spdLength<<","<<spdWidth<<","<<fSegmentation->Dy()<<","<<fNPixelsX<<","<<fNPixelsZ<<","<<xPitch<<","<<zPitch<<endl;
+    cout<<"len,wid,thickness,nx,nz,pitchx,pitchz,difcoef ="<<spdLength<<","<<spdWidth<<","<<spdThickness<<","<<fNPixelsX<<","<<fNPixelsZ<<","<<xPitch<<","<<zPitch<<","<<difCoef<<endl;
   //  Array of pointers to the label-signal list
 
     Int_t maxNDigits = fNPixelsX*fNPixelsZ + fNPixelsX ;; 
@@ -143,21 +144,26 @@ void AliITSsimulationSPDdubna::DigitiseModule(AliITSmodule *mod, Int_t module, I
     static Bool_t first;
     Int_t lasttrack=-2;
     Int_t hit, iZi, jz, jx;
-    //cout<<"SPD: module,nhits ="<<module<<","<<nhits<<endl;
-    Int_t idhit=-1;
+    Int_t idhit=-1; //!
+    cout<<"SPDdubna: module,nhits ="<<module<<","<<nhits<<endl;
     for (hit=0;hit<nhits;hit++) {
         AliITShit *iHit = (AliITShit*) fHits->At(hit);
-       Int_t layer = iHit->GetLayer();
-        Float_t yPix0 = -73; 
-        if(layer == 1) yPix0 = -77; 
-
+       //Int_t layer = iHit->GetLayer();
+        Float_t yPix0 = -spdThickness/2; 
+
+       // work with the idtrack=entry number in the TreeH
+       //Int_t idhit,idtrack; //!
+       //mod->GetHitTrackAndHitIndex(hit,idtrack,idhit);  //!    
+       //Int_t idtrack=mod->GetHitTrackIndex(hit);  
+        // or store straight away the particle position in the array
+       // of particles : 
        if(iHit->StatusEntering()) idhit=hit;
         Int_t itrack = iHit->GetTrack();
         Int_t dray = 0;
    
        if (lasttrack != itrack || hit==(nhits-1)) first = kTRUE; 
 
-       //        Int_t parent = iHit->GetParticle()->GetFirstMother();
+       //Int_t parent = iHit->GetParticle()->GetFirstMother();
         Int_t partcode = iHit->GetParticle()->GetPdgCode();
 
 //  partcode (pdgCode): 11 - e-, 13 - mu-, 22 - gamma, 111 - pi0, 211 - pi+
@@ -188,28 +194,23 @@ void AliITSsimulationSPDdubna::DigitiseModule(AliITSmodule *mod, Int_t module, I
 
        // Get track status
        Int_t status = iHit->GetTrackStatus();      
-       //cout<<"hit,status,y ="<<hit<<","<<status<<","<<yPix<<endl;      
 
        // Check boundaries
        if(zPix  > spdLength/2) {
-         //cout<<"!!!1 z outside ="<<zPix<<endl;
+         //cout<<"!!! SPD: z outside ="<<zPix<<endl;
          zPix = spdLength/2 - 10;
-        //cout<<"!!!2 z outside ="<<zPix<<endl;
        }
        if(zPix  < 0 && zPix < -spdLength/2) {
-         //cout<<"!!!1 z outside ="<<zPix<<endl;
+         //cout<<"!!! SPD: z outside ="<<zPix<<endl;
          zPix = -spdLength/2 + 10;
-        //cout<<"!!!2 z outside ="<<zPix<<endl;
        }
        if(xPix  > spdWidth/2) {
-         //cout<<"!!!1 x outside ="<<xPix<<endl;
+         //cout<<"!!! SPD: x outside ="<<xPix<<endl;
          xPix = spdWidth/2 - 10;
-        //cout<<"!!!2 x outside ="<<xPix<<endl;
        }
        if(xPix  < 0 && xPix < -spdWidth/2) {
-         //cout<<"!!!1 x outside ="<<xPix<<endl;
+         //cout<<"!!! SPD: x outside ="<<xPix<<endl;
          xPix = -spdWidth/2 + 10;
-        //cout<<"!!!2 x outside ="<<xPix<<endl;
        }
        Int_t trdown = 0;
 
@@ -259,7 +260,7 @@ void AliITSsimulationSPDdubna::DigitiseModule(AliITSmodule *mod, Int_t module, I
        // number of the steps along the track:
        Int_t nsteps = ndZ;
        if(ndX > ndZ) nsteps = ndX;
-       if(nsteps < 6) nsteps = 6;  // minimum number of the steps 
+       if(nsteps < 20) nsteps = 20;  // minimum number of the steps 
 
        if (projDif < 5 ) {
           drPath = (yPix-yPix0)*1.e-4;  
@@ -523,7 +524,7 @@ void AliITSsimulationSPDdubna::ChargeToSignal(Float_t **pList)
 
   AliITS *aliITS = (AliITS*)gAlice->GetModule("ITS");
   
-
   Float_t threshold = (float)fResponse->MinVal();
 
   Int_t digits[3], tracks[3], hits[3],gi,j1;
@@ -591,8 +592,8 @@ void AliITSsimulationSPDdubna::CreateHistograms()
       printf("SPD - create histograms\n");
 
       fHis=new TObjArray(fNPixelsZ);
+      TString spdName("spd_");
       for (Int_t i=0;i<fNPixelsZ;i++) {
-       TString spdName("spd_");
           Char_t pixelz[4];
           sprintf(pixelz,"%d",i+1);
           spdName.Append(pixelz);
index fefcd02..d7a5df9 100644 (file)
@@ -183,7 +183,24 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
      printf("Found %d entries in the tree (must be one per module per event!)\n",nent);
      Int_t lay, lad, det;
      AliITSgeom *geom = ITS->GetITSgeom();
-   
+
+   AliITSDetType *iDetType=ITS->DetType(0);
+   AliITSsegmentationSPD *seg0=(AliITSsegmentationSPD*)iDetType->GetSegmentationModel();
+   printf("SPD dimensions %f %f %f \n",seg0->Dx(),seg0->Dz(),seg0->Dy());
+   printf("SPD npixels %d %d \n",seg0->Npz(),seg0->Npx());
+   printf("SPD pitches %d %d \n",seg0->Dpz(0),seg0->Dpx(0));
+   Float_t SPDlength = seg0->Dz();     
+   Float_t SPDwidth = seg0->Dx();      
+   Float_t SPDthickness = seg0->Dy();
+   Float_t xpitch = seg0->Dpx(0);
+   if(SPDlength > 80000) SPDthickness = 150;
+   Float_t ylim;
+   if(SPDthickness < 200) {
+     ylim = SPDthickness/2 - 4;
+   }else{
+     ylim = SPDthickness/2 - 10;
+   }
+
      for (Int_t idettype=0;idettype<3;idettype++) {
 
        TClonesArray *ITSclusters  = ITS->ClustersAddress(idettype);
@@ -218,18 +235,26 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
 
                Int_t clustersizex = itsclu->NclX();
                Int_t clustersizez = itsclu->NclZ();
-               //       Int_t xstart = itsclu->XStart();
-               //       Int_t xstop = itsclu->XStop();
                Int_t xstart = itsclu->XStartf();
                Int_t xstop = itsclu->XStopf();
-               Float_t fxstart = xstart*50;
-               Float_t fxstop = (xstop+1)*50;
+               Float_t fxstart = xstart*xpitch;
+               Float_t fxstop = (xstop+1)*xpitch;
                Float_t zstart = itsclu->ZStart();
                Float_t zstop = itsclu->ZStop();
                Int_t zend = itsclu->Zend();
                Int_t ntrover = itsclu->NTracks();
-               Float_t clusterx = itsclu->X();
+                Float_t clusterx = 0;
                Float_t clusterz = itsclu->Z();
+
+               Int_t tr0, tr1, tr2;
+               itsclu->GetTracks(tr0,tr1,tr2);
+
+               for(Int_t ii=0;ii<clustersizex;ii++) {
+                 clusterx += (xstart+0.5+ii)*xpitch;
+               }
+               clusterx /= clustersizex;
+               clusterz /= clustersizez;   
+
                Float_t clusterQ = itsclu->Q();
 
                if(lay == 1) occup1 += clusterQ;                
@@ -248,11 +273,8 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
                Float_t dxprimlast = 10.e+6;
                Float_t dzprimlast = 10.e+6;
 
-               Float_t SPDlength = 83600;      
-               Float_t SPDwidth = 12800;       
-                Float_t xhit0 = 1e+5;
-                Float_t zhit0 = 1e+5;
-               
+                Float_t xhit0 = 1e+6;
+                Float_t zhit0 = 1e+6;
        for (Int_t hit=0;hit<nhits;hit++) {
 
                  // Find coordinate differences between the hit and cluster positions
@@ -267,9 +289,16 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
                  Int_t track = itsHit->GetTrack();
                  Int_t dray = 0;
                  Int_t hitstat = itsHit->GetTrackStatus();
-
                  Float_t zhit = 10000*itsHit->GetZL();
                  Float_t xhit = 10000*itsHit->GetXL();
+                 Float_t yhit = 10000*itsHit->GetYL();
+        Int_t parent = itsHit->GetParticle()->GetFirstMother();
+        Int_t partcode = itsHit->GetParticle()->GetPdgCode();
+
+        Int_t hitprim = 0;
+
+        if(parent < 0) hitprim = hitprim + 1; // hitprim=1 for the primery
+       // particles
 
         Float_t pxsimL = itsHit->GetPXL();  // the momenta at GEANT points
         Float_t pysimL = itsHit->GetPYL();
@@ -296,41 +325,49 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
 
                  zhit += SPDlength/2;
                  xhit += SPDwidth/2;
-                 Float_t yhit = 10000*itsHit->GetYL();
 
-                 if(hitlayer == 1 && hitstat == 66 && yhit > 71) {
+
+
+                 if(hitlayer == 1 && hitstat == 66 && yhit > ylim) {
+                   //if(hitstat == 66 && hitprim ==1) {
                    xhit0 = xhit;
                    zhit0 = zhit;
                  }
-                 if(hitlayer == 2 && hitstat == 66 && yhit < -71) {
+
+                 
+                 if(hitlayer == 2 && hitstat == 66 && yhit < -ylim) {
                    xhit0 = xhit;
                    zhit0 = zhit;
-                 }
-
+                 }              
+                         
                  if(hitstat != 68) continue; // Take only the hit if the last
                  // track point went out from
                  // the detector.
 
-                 if(xhit0 > 9e+4 || zhit0 > 9e+4) continue;
+                 if(xhit0 > 1e+5 || zhit0 > 1e+5) continue;
 
                  Float_t xmed = (xhit + xhit0)/2;
                  Float_t zmed = (zhit + zhit0)/2;
 
                  Float_t xdif = xmed - clusterx;
                  Float_t zdif = zmed - clusterz;
+                 xhit0 = 1e+6;
+                 zhit0 = 1e+6;
+
+        // Consider the hits inside of cluster region only   b.b.
 
+                 if((xmed >= fxstart && xmed <= fxstop) && (zmed >= zstart && zmed <= zstop)) {
 
-        // Consider the hits inside of cluster region only
+         // Consider the hits only with the track number equaled to one
+         // of the cluster
+                   //if((track == tr0) || (track == tr1) || (track == tr2)) {      
 
-  if((xmed >= fxstart && xmed <= fxstop) && (zmed >= zstart && zmed <= zstop)) {
         icl = 1;
 
                 //        part = (TParticle *)particles.UncheckedAt(track);
                 //        Int_t partcode = part->GetPdgCode();
                 //              Int_t primery = gAlice->GetPrimary(track);
 
-        Int_t parent = itsHit->GetParticle()->GetFirstMother();
-        Int_t partcode = itsHit->GetParticle()->GetPdgCode();
 
 //  partcode (pdgCode): 11 - e-, 13 - mu-, 22 - gamma, 111 - pi0, 211 - pi+
 //                      310 - K0s, 321 - K+, 2112 - n, 2212 - p, 3122 - lambda
@@ -352,9 +389,10 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
        Float_t theta = itsHit->GetParticle()->Theta(); // Theta angle at the
                                                   // vertex
        //Float_t eta = itsHit->GetParticle()->Eta(); // Pseudo rapidity at the
-                                                  // vertex
+                                                          // vertex
+       Float_t y;
        if((energy-pz) > 0) {
-         Float_t y = 0.5*TMath::Log((energy+pz)/(energy-pz));
+         y = 0.5*TMath::Log((energy+pz)/(energy-pz));
        }else{
          cout<<" Warning: energy < pz ="<<energy<<","<<pz<<endl;
          y = 10;
@@ -368,7 +406,7 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
         Float_t pzsim = itsHit->GetPZG();
        Float_t psim = TMath::Sqrt(pxsim*pxsim+pysim*pysim+pzsim*pzsim);
 
-        Int_t hitprim = 0;
+        //Int_t hitprim = 0;
 
         if(partcode == 11 && pmod < 6) dray = 1; // delta ray is e-
                                                  // at p < 6 MeV/c
@@ -379,7 +417,7 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
                                                  // detector and returned
                                                  // again
 
-        if(parent < 0) hitprim = hitprim + 1; // hitprim=1 for the primery
+        //if(parent < 0) hitprim = hitprim + 1; // hitprim=1 for the primery
                                               // particles
 
         if(hitprim > 0) noverprim = noverprim + 1;
@@ -443,15 +481,13 @@ void SPDclusterTestDubna (Int_t evNumber1=0,Int_t evNumber2=0)
       } // primery particles
 
      } // end of cluster region
+
    } // end of hit loop
 
       if(icl == 1) {
 
         // fill ntuple1
 
-        //      ntuple1->Fill(clusterlayer,clustersizex,clustersizez,noverlaps,\
-noverprim,dx,dz);
-
         if(noverlaps == 0) noverlaps = 1; // cluster contains one or more
                                           // delta rays only