]> git.uio.no Git - u/mrichter/AliRoot.git/blobdiff - ZDC/AliZDCv1.cxx
Correction for superposition of ZDC volumes with MUON arm one
[u/mrichter/AliRoot.git] / ZDC / AliZDCv1.cxx
index 820a4742c2ca1b1b3a1f19d9db329816ace8e4ec..b1c3ea6debe1a530ac6141be7edb0b9cda9b9c44 100644 (file)
 
 /*
 $Log$
+Revision 1.16  2001/03/15 16:12:04  coppedis
+Code review
+
+Revision 1.15  2001/03/12 17:47:56  hristov
+Changes needed on Sun with CC 5.0
+
+Revision 1.14  2001/02/23 16:48:28  coppedis
+Correct bug in ZEM hit definition
+
+Revision 1.13  2001/02/07 18:07:41  coppedis
+Modif for splitting
+
+Revision 1.12  2001/01/26 19:56:27  hristov
+Major upgrade of AliRoot code
+
+Revision 1.11  2001/01/16 07:43:33  hristov
+Initialisation of ZDC hits
+
+Revision 1.10  2000/12/14 15:20:02  coppedis
+Hits2Digits method for digitization
+
 Revision 1.9  2000/12/13 10:33:49  coppedis
 Prints only if fDebug==1
 
@@ -43,7 +64,7 @@ Revision 1.2  2000/07/11 11:12:34  fca
 Some syntax corrections for non standard HP aCC
 
 Revision 1.1  2000/07/10 13:58:01  fca
-New version of ZDC from E.Scomparin & C.Oppedisano
+New version of ZDC from E.Scomparin & C.Oppedisano
 
 Revision 1.7  2000/01/19 17:17:40  fca
 
@@ -107,10 +128,9 @@ AliZDCv1::AliZDCv1() : AliZDC()
   fMedSensF2  = 0;
   fMedSensZN  = 0;
   fMedSensZP  = 0;
-  fMedSensGR  = 0;
   fMedSensZEM = 0;
+  fMedSensGR  = 0;
   fMedSensPI  = 0;
-  fNoShower   = 0;
 }
  
 //_____________________________________________________________________________
@@ -121,16 +141,61 @@ AliZDCv1::AliZDCv1(const char *name, const char *title)
   // Standard constructor for Zero Degree Calorimeter 
   //
 
-  fDigits = new TClonesArray("AliZDCDigit",1000);
-
   fMedSensF1  = 0;
   fMedSensF2  = 0;
   fMedSensZN  = 0;
   fMedSensZP  = 0;
-  fMedSensGR  = 0;
   fMedSensZEM = 0;
+  fMedSensGR  = 0;
   fMedSensPI  = 0;
-  fNoShower   = 0;
+
+  
+  // Parameters for light tables
+  fNalfan = 90;       // Number of Alfa (neutrons)
+  fNalfap = 90;       // Number of Alfa (protons)
+  fNben = 18;         // Number of beta (neutrons)
+  fNbep = 28;         // Number of beta (protons)
+  Int_t ip,jp,kp;
+  for(ip=0; ip<4; ip++){
+     for(kp=0; kp<fNalfap; kp++){
+        for(jp=0; jp<fNbep; jp++){
+           fTablep[ip][kp][jp] = 0;
+        } 
+     }
+  }
+  Int_t in,jn,kn;
+  for(in=0; in<4; in++){
+     for(kn=0; kn<fNalfan; kn++){
+        for(jn=0; jn<fNben; jn++){
+           fTablen[in][kn][jn] = 0;
+        } 
+     }
+  }
+
+  // Parameters for hadronic calorimeters geometry
+  fDimZP[0] = 11.2;
+  fDimZP[1] = 6.;
+  fDimZP[2] = 75.;    
+  fPosZN[0] = 0.;
+  fPosZN[1] = 0.;
+  fPosZN[2] = 11650.;
+  fPosZP[0] = -23.;
+  fPosZP[1] = 0.;
+  fPosZP[2] = 11600.;
+  fFibZN[0] = 0.;
+  fFibZN[1] = 0.01825;
+  fFibZN[2] = 50.;
+  fFibZP[0] = 0.;
+  fFibZP[1] = 0.0275;
+  fFibZP[2] = 75.;
+  
+  // Parameters for EM calorimeter geometry
+  fPosZEM[0] = 0.;
+  fPosZEM[1] = 5.8;
+  fPosZEM[2] = 11600.;
+  
+
+  fDigits = new TClonesArray("AliZDCDigit",1000);
 }
  
 //_____________________________________________________________________________
@@ -149,20 +214,17 @@ void AliZDCv1::CreateGeometry()
 void AliZDCv1::CreateBeamLine()
 {
   
-  Float_t angle;
-  Float_t zq, conpar[9], elpar[3], tubpar[3];
+  Float_t angle, zq, conpar[9], elpar[3], tubpar[3], zd1, zd2;
   Int_t im1, im2;
-  Float_t zd1, zd2;
-  
   
   Int_t *idtmed = fIdtmed->GetArray();
   
-  // -- Mother of the ZDC 
+  // -- Mother of the ZDCs 
   
   conpar[0] = 0.;
   conpar[1] = 360.;
   conpar[2] = 2.;
-  conpar[3] = 805.;
+  conpar[3] = 2000.;
   conpar[4] = 0.;
   conpar[5] = 55.;
   conpar[6] = 13060.;
@@ -174,20 +236,21 @@ void AliZDCv1::CreateBeamLine()
   // -- FIRST SECTION OF THE BEAM PIPE (from compensator dipole to 
   //    beginning of D1) 
   
-  zd1 = 1921.6;
+  zd1 = 2000.;
   
   tubpar[0] = 6.3/2.;
   tubpar[1] = 6.7/2.;
-  tubpar[2] = 3916.7/2.;
+  tubpar[2] = 3838.3/2.;
   gMC->Gsvolu("P001", "TUBE", idtmed[5], tubpar, 3);
   gMC->Gspos("P001", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
   
   //-- SECOND SECTION OF THE BEAM PIPE (FROM THE END OF D1 TO THE BEGINNING OF
   //    D2) 
   
-  //-- FROM MAGNETIC BEGINNING OG D1 TO MAGNETIC END OF D1 + 23.5 cm
-  //-- Elliptic pipe
+  //-- FROM MAGNETIC BEGINNING OF D1 TO MAGNETIC END OF D1 + 23.5 cm
+  //-- Elliptic pipe -> ELLIPTIC PIPE NOT INSERTED FOR THE MOMENT!
   
+  //-> D1 begins at (6310.8-472.5)
   zd1 = 6310.8-472.5;
   
   elpar[0] = 6.84/2.;
@@ -234,7 +297,7 @@ void AliZDCv1::CreateBeamLine()
   gMC->Gsvolu("P002", "TUBE", idtmed[5], tubpar, 3);
   gMC->Gspos("P002", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
   
-  zd1 += tubpar[2] * 2.;
+  zd1 += tubpar[2]*2.;
   
   tubpar[0] = 10./2.;
   tubpar[1] = 10.4/2.;
@@ -405,6 +468,9 @@ void AliZDCv1::CreateBeamLine()
   // --  END OF BEAM PIPE VOLUME DEFINITION. MAGNET DEFINITION FOLLOWS 
   //     (LHC OPTICS 6) 
   
+  // ----------------------------------------------------------------
+  //                   Replaced by the muon dipole
+  // ----------------------------------------------------------------
   // -- COMPENSATOR DIPOLE (MBXW) 
   //     GAP (VACUUM WITH MAGNETIC FIELD) 
   
@@ -422,22 +488,25 @@ void AliZDCv1::CreateBeamLine()
 //  gMC->Gsvolu("YMBX", "TUBE", idtmed[5], tubpar, 3);
 //  gMC->Gspos("YMBX", 1, "ZDC ", 0., 0., tubpar[2] + 805., 0, "ONLY");
   
+  // ----------------------------------------------------------------
+  //                 Replaced by the muon compesator
+  // ----------------------------------------------------------------
   // -- COMPENSATOR DIPOLE (MCBWA) 
   //     GAP (VACUUM WITH MAGNETIC FIELD) 
   
-  tubpar[0] = 0.;
-  tubpar[1] = 4.5;
-  tubpar[2] = 170./2.;
-  gMC->Gsvolu("MCBW", "TUBE", idtmed[11], tubpar, 3);
-  gMC->Gspos("MCBW", 1, "ZDC ", 0., 0., tubpar[2] + 1921.6, 0, "ONLY");
+//  tubpar[0] = 0.;
+//  tubpar[1] = 4.5;
+//  tubpar[2] = 170./2.;
+//  gMC->Gsvolu("MCBW", "TUBE", idtmed[11], tubpar, 3);
+//  gMC->Gspos("MCBW", 1, "ZDC ", 0., 0., tubpar[2] + 1921.6, 0, "ONLY");
   
   // --  YOKE (IRON WITHOUT MAGNETIC FIELD) 
   
-  tubpar[0] = 4.5;
-  tubpar[1] = 55.;
-  tubpar[2] = 170./2.;
-  gMC->Gsvolu("YMCB", "TUBE", idtmed[5], tubpar, 3);
-  gMC->Gspos("YMCB", 1, "ZDC ", 0., 0., tubpar[2] + 1921.6, 0, "ONLY");
+//  tubpar[0] = 4.5;
+//  tubpar[1] = 55.;
+//  tubpar[2] = 170./2.;
+//  gMC->Gsvolu("YMCB", "TUBE", idtmed[5], tubpar, 3);
+//  gMC->Gspos("YMCB", 1, "ZDC ", 0., 0., tubpar[2] + 1921.6, 0, "ONLY");
   
   // -- INNER TRIPLET 
   
@@ -538,9 +607,31 @@ void AliZDCv1::CreateBeamLine()
 void AliZDCv1::CreateZDC()
 {
   
-  Int_t *idtmed = fIdtmed->GetArray();
   Int_t irot1, irot2;
   Float_t DimPb[6], DimVoid[6];
+  
+  Int_t *idtmed = fIdtmed->GetArray();
+
+  // Parameters for hadronic calorimeters geometry
+  // NB -> parameters used ONLY in CreateZDC()
+  Float_t fDimZN[3] = {3.52, 3.52, 50.};  // Dimensions of neutron detector
+  Float_t fGrvZN[3] = {0.03, 0.03, 50.};  // Grooves for neutron detector
+  Float_t fGrvZP[3] = {0.04, 0.04, 75.};  // Grooves for proton detector
+  Int_t   fDivZN[3] = {11, 11, 0};       // Division for neutron detector
+  Int_t   fDivZP[3] = {7, 15, 0};        // Division for proton detector
+  Int_t   fTowZN[2] = {2, 2};                    // Tower for neutron detector
+  Int_t   fTowZP[2] = {4, 1};                    // Tower for proton detector
+
+  // Parameters for EM calorimeter geometry
+  // NB -> parameters used ONLY in CreateZDC()
+  Float_t fDimZEMPb  = 0.15*(TMath::Sqrt(2.));  // z-dimension of the Pb slice
+  Float_t fDimZEMAir = 0.001;                  // scotch
+  Float_t fFibRadZEM = 0.0315;                         // External fiber radius (including cladding)
+  Int_t   fDivZEM[3] = {92, 0, 20};            // Divisions for EM detector
+  Float_t fDimZEM0 = 2*fDivZEM[2]*(fDimZEMPb+fDimZEMAir+fFibRadZEM*(TMath::Sqrt(2.)));
+  Float_t fDimZEM[6] = {fDimZEM0, 3.5, 3.5, 45., 0., 0.}; // Dimensions of EM detector
+  Float_t fFibZEM2 = fDimZEM[2]/TMath::Sin(fDimZEM[3]*kDegrad)-fFibRadZEM;
+  Float_t fFibZEM[3] = {0., 0.0275, fFibZEM2};  // Fibers for EM calorimeter
 
   
   //-- Create calorimeters geometry
@@ -812,12 +903,11 @@ void AliZDCv1::CreateMaterials()
   Int_t *idtmed = fIdtmed->GetArray();
   
   Float_t dens, ubuf[1], wmat[2], a[2], z[2], epsil=0.001, stmin=0.01;
-  Int_t   i, isvolActive, isvol, inofld;
+  Float_t deemax = -1, stemax;
   Float_t fieldm = gAlice->Field()->Max();
   Float_t tmaxfd=gAlice->Field()->Max();
+  Int_t   i, isvolActive, isvol, inofld;
   Int_t   isxfld = gAlice->Field()->Integ();
-  Float_t deemax=-1;
-  Float_t stemax;
   
   // --- Store in UBUF r0 for nuclear radius calculation R=r0*A**1/3 
 
@@ -1001,12 +1091,7 @@ void AliZDCv1::Init()
 void AliZDCv1::InitTables()
 {
   Int_t k, j;
-  //Initialize parameters for light tables and read them
-  fNalfan = 90;
-  fNalfap = 90;
-  fNben = 18;
-  fNbep = 28;
-  
+
   char *lightfName1,*lightfName2,*lightfName3,*lightfName4,
        *lightfName5,*lightfName6,*lightfName7,*lightfName8;
   FILE *fp1, *fp2, *fp3, *fp4, *fp5, *fp6, *fp7, *fp8;
@@ -1089,39 +1174,33 @@ Int_t AliZDCv1::Digitize(Int_t Det, Int_t Quad, Int_t Light)
     printf("\n Digitize -> Det = %d, Quad = %d, Light = %d\n", Det, Quad, Light);
   }   
   
+  // Parameters for conversion of light yield in ADC channels
+  Float_t fPMGain[3][5];      // PM gain
+  Float_t fADCRes;            // ADC conversion factor
+  
   Int_t j,i;
   for(i=0; i<3; i++){
      for(j=0; j<5; j++){
-        fPedMean[i][j]  = 50.;
-        fPedSigma[i][j] = 10.;
-        fPMGain[i][j]   = 10000000.;
+        fPMGain[i][j]   = 100000.;
      }
   }
   fADCRes   = 0.00000064; // ADC Resolution: 250 fC/ADCch
   
-  Float_t Ped = gRandom->Gaus(fPedMean[Det-1][Quad],fPedSigma[Det-1][Quad]);
-  Int_t ADCch = Int_t(Light*fPMGain[Det-1][Quad]*fADCRes+Ped);
-  
-  if(fDebug == 1){
-    printf("   Ped = %f, ADCch = %d\n", Ped, ADCch);
-  }  
-   
+  Int_t ADCch = Int_t(Light*fPMGain[Det-1][Quad]*fADCRes);
+     
   return ADCch;
 }
+
+
 //_____________________________________________________________________________
-void AliZDCv1::Hits2Digits(Int_t ntracks)
+void AliZDCv1::SDigits2Digits()
 {
-  // Creation of the digits from hits 
-
-  if(fDigits!=0) fDigits->Clear();
-  else fDigits = new TClonesArray ("AliZDCDigit",1000);
+   Hits2Digits(gAlice->GetNtrack());
+}
 
-  char branchname[10];
-  sprintf(branchname,"%s",GetName());
-  gAlice->TreeD()->Branch(branchname,&fDigits, fBufferSize);
-  
-  gAlice->TreeD()->GetEvent(0);
-  
+//_____________________________________________________________________________
+void AliZDCv1::Hits2Digits(Int_t ntracks)
+{
   AliZDCDigit *newdigit;
   AliZDCHit   *hit;
 
@@ -1162,6 +1241,8 @@ void AliZDCv1::Hits2Digits(Int_t ntracks)
         }
      } // Hits loop
   
+  } // Tracks loop
+  
      if(fDebug == 1){
        printf("\n       PMCZN = %d, PMQZN[0] = %d, PMQZN[1] = %d, PMQZN[2] = %d, PMQZN[3] = %d\n"
            , PMCZN, PMQZN[0], PMQZN[1], PMQZN[2], PMQZN[3]);
@@ -1204,8 +1285,6 @@ void AliZDCv1::Hits2Digits(Int_t ntracks)
      new((*fDigits)[fNdigits]) AliZDCDigit(*newdigit);
      fNdigits++;
      delete newdigit;
-  
-  } // Tracks loop
       
   
   gAlice->TreeD()->Fill();
@@ -1218,7 +1297,7 @@ void AliZDCv1::Hits2Digits(Int_t ntracks)
   
 }
 //_____________________________________________________________________________
- void AliZDCv1::MakeBranch(Option_t *opt)
+ void AliZDCv1::MakeBranch(Option_t *opt, char *file)
 {
   //
   // Create a new branch in the current Root Tree
@@ -1226,12 +1305,20 @@ void AliZDCv1::Hits2Digits(Int_t ntracks)
 
   AliDetector::MakeBranch(opt);
   
-  char branchname[10];
+  Char_t branchname[10];
   sprintf(branchname,"%s",GetName());
-  char *cD = strstr(opt,"D");
+  const char *cD = strstr(opt,"D");
+
+  if (gAlice->TreeD() && cD) {
+
+    // Creation of the digits from hits 
 
-  if (fDigits   && gAlice->TreeD() && cD) {
-    gAlice->TreeD()->Branch(branchname,&fDigits, fBufferSize);
+    if(fDigits!=0) fDigits->Clear();
+    else fDigits = new TClonesArray ("AliZDCDigit",1000);
+    char branchname[10];
+    sprintf(branchname,"%s",GetName());
+    gAlice->MakeBranchInTree(gAlice->TreeD(), 
+                             branchname, &fDigits, fBufferSize, file) ;
     printf("* AliZDCv1::MakeBranch    * Making Branch %s for digits\n\n",branchname);
   }     
 }
@@ -1243,12 +1330,12 @@ void AliZDCv1::StepManager()
   //
 
   Int_t j;
-
   Int_t vol[2], ibeta=0, ialfa, ibe, nphe;
   Float_t x[3], xdet[3], destep, hits[10], m, ekin, um[3], ud[3], be, radius, out;
   TLorentzVector s, p;
   const char *knamed;
 
+  for (j=0;j<10;j++) hits[j]=0;
 
   if((gMC->GetMedium() == fMedSensZN) || (gMC->GetMedium() == fMedSensZP) ||
      (gMC->GetMedium() == fMedSensGR) || (gMC->GetMedium() == fMedSensF1) ||
@@ -1284,10 +1371,10 @@ void AliZDCv1::StepManager()
     if(vol[0]==1){
       xdet[0] = x[0]-fPosZN[0];
       xdet[1] = x[1]-fPosZN[1];
-      if((xdet[0]<=0.) && (xdet[1]>=0.)) vol[1]=1;
-      if((xdet[0]>0.) && (xdet[1]>0.)) vol[1]=2;
-      if((xdet[0]<0.) && (xdet[1]<0.)) vol[1]=3;
-      if((xdet[0]>0.) && (xdet[1]<0.)) vol[1]=4;
+      if((xdet[0]<=0.) && (xdet[1]>=0.))  vol[1]=1;
+      if((xdet[0]>0.)  && (xdet[1]>0.))   vol[1]=2;
+      if((xdet[0]<0.)  && (xdet[1]<0.))   vol[1]=3;
+      if((xdet[0]>0.)  && (xdet[1]<0.))   vol[1]=4;
     }
     
     //Quadrant in ZP
@@ -1301,8 +1388,8 @@ void AliZDCv1::StepManager()
          if(xqZP>=(i-3) && xqZP<(i-2)){
           vol[1] = i;
           break;
-       }
-     }
+        }
+      }
     }
     
     //ZEM has only 1 quadrant
@@ -1310,25 +1397,15 @@ void AliZDCv1::StepManager()
       vol[1] = 1;
       xdet[0] = x[0]-fPosZEM[0];
       xdet[1] = x[1]-fPosZEM[1];
-//      printf("x %f %f xdet %f %f\n",x[0],x[1],xdet[0],xdet[1]);
     }
 
-//    if(vol[1]>4){
-//    printf("\n-> Det. %d Quad. %d \n", vol[0], vol[1]);
-//    printf("x %f %f xdet %f %f\n",x[0],x[1],xdet[0],xdet[1]);}
 
   // Store impact point and kinetic energy of the ENTERING particle
     
-//    Int_t Curtrack = gAlice->CurrentTrack();
-//    Int_t Prim = gAlice->GetPrimary(Curtrack);
-//    printf ("Primary: %d, Current Track: %d \n", Prim, Curtrack); 
-    
 //    if(Curtrack==Prim){
       if(gMC->IsTrackEntering()){
         //Particle energy
         gMC->TrackMomentum(p);
-//      printf("p[0] = %f, p[1] = %f, p[2] = %f, p[3] = %f \n", 
-//                 p[0], p[1], p[2], p[3]);
         hits[3] = p[3];
 
         // Impact point on ZDC  
@@ -1385,8 +1462,6 @@ void AliZDCv1::StepManager()
        gMC->TrackMomentum(p);
        Float_t ptot=TMath::Sqrt(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);
        Float_t beta =  ptot/p[3];
-//       Int_t pcID = gMC->TrackPid();
-//       printf("      Pc %d in quadrant %d -> beta = %f \n", pcID, vol[1], beta);
        if(beta<0.67) return;
        if((beta>=0.67) && (beta<=0.75)) ibeta = 0;
        if((beta>0.75)  && (beta<=0.85)) ibeta = 1;
@@ -1469,8 +1544,8 @@ void AliZDCv1::StepManager()
          if(ibe>fNbep) ibe=fNbep;
          out =  charge*charge*fTablep[ibeta][ialfa][ibe];
         nphe = gRandom->Poisson(out);
-        hits[7] = nphe;        //fLightPMQ
-        hits[8] = 0;
+        hits[7] = 0;   
+        hits[8] = nphe;        //fLightPMC
         hits[9] = 0;
         AddHit(gAlice->CurrentTrack(), vol, hits);
        }