ZDC/AliZDCv1.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /*
  17 $Log$
  18 Revision 1.4  2000/08/24 09:23:59  hristov
  19 Bug in ZDC geometry corrected by E.Scomparin
  20
  21 Revision 1.3  2000/07/12 06:59:16  fca
  22 Fixing dimension of hits array
  23
  24 Revision 1.2  2000/07/11 11:12:34  fca
  25 Some syntax corrections for non standard HP aCC
  26
  27 Revision 1.1  2000/07/10 13:58:01  fca
  28 New version of ZDC from E.Scomparin & C.Oppedisano
  29
  30 Revision 1.7  2000/01/19 17:17:40  fca
  31
  32 Revision 1.6  1999/09/29 09:24:35  fca
  33 Introduction of the Copyright and cvs Log
  34
  35 */
  36
  37 ///////////////////////////////////////////////////////////////////////////////
  38 //                                                                           //
  39 //  Zero Degree Calorimeter                                                  //
  40 //  This class contains the basic functions for the Time Of Flight           //
  41 //  detector. Functions specific to one particular geometry are              //
  42 //  contained in the derived classes                                         //
  43 //                                                                           //
  44 ///////////////////////////////////////////////////////////////////////////////
  45
  46 #include <TBRIK.h>
  47 #include <TNode.h>
  48 #include <TMath.h>
  49
  50 #include "stdio.h"
  51 #include "AliZDCv1.h"
  52 #include "AliRun.h"
  53 #include "AliMagF.h"
  54 #include "AliMC.h"
  55 #include "AliCallf77.h"
  56 #include "AliConst.h"
  57 #include "AliPDG.h"
  58
  59
  60 ClassImp(AliZDCv1)
  61
  62
  63 ///////////////////////////////////////////////////////////////////////////////
  64 //                                                                           //
  65 //  Zero Degree Calorimeter version 1                                        //
  66 //                                                                           //
  67 ///////////////////////////////////////////////////////////////////////////////
  68
  69 //_____________________________________________________________________________
  70 AliZDCv1::AliZDCv1() : AliZDC()
  71 {
  72   //
  73   // Default constructor for Zero Degree Calorimeter
  74   //
  75   fMedSensF1 = 0;
  76   fMedSensF2 = 0;
  77   fMedSensZN = 0;
  78   fMedSensZP = 0;
  79   fMedSensGR = 0;
  80 }
  81
  82 //_____________________________________________________________________________
  83 AliZDCv1::AliZDCv1(const char *name, const char *title)
  84   : AliZDC(name,title)
  85 {
  86   //
  87   // Standard constructor for Zero Degree Calorimeter
  88   //
  89   fMedSensF1 = 0;
  90   fMedSensF2 = 0;
  91   fMedSensZN = 0;
  92   fMedSensZP = 0;
  93   fMedSensGR = 0;
  94 }
  95
  96 //_____________________________________________________________________________
  97 void AliZDCv1::CreateGeometry()
  98 {
  99   //
 100   // Create the geometry for the Zero Degree Calorimeter version 1
 101   //* Initialize COMMON block ZDC_CGEOM
 102   //*
 103
 104   CreateBeamLine();
 105   CreateZDC();
 106 }
 107
 108 //_____________________________________________________________________________
 109 void AliZDCv1::CreateBeamLine()
 110 {
 111
 112   Float_t angle;
 113   Float_t zq, conpar[9], elpar[3], tubpar[3];
 114   Int_t im1, im2;
 115   Float_t zd1, zd2;
 116
 117
 118   Int_t *idtmed = fIdtmed->GetArray();
 119
 120   // -- Mother of the ZDC
 121
 122   conpar[0] = 0.;
 123   conpar[1] = 360.;
 124   conpar[2] = 2.;
 125   conpar[3] = 805.;
 126   conpar[4] = 0.;
 127   conpar[5] = 55.;
 128   conpar[6] = 13060.;
 129   conpar[7] = 0.;
 130   conpar[8] = 55.;
 131   gMC->Gsvolu("ZDC ", "PCON", idtmed[10], conpar, 9);
 132   gMC->Gspos("ZDC ", 1, "ALIC", 0., 0., 0., 0, "ONLY");
 133
 134   // -- FIRST SECTION OF THE BEAM PIPE (from compensator dipole to
 135   //    beginning of D1)
 136
 137   zd1 = 1921.6;
 138
 139   tubpar[0] = 6.3/2.;
 140   tubpar[1] = 6.7/2.;
 141   tubpar[2] = 3916.7/2.;
 142   gMC->Gsvolu("P001", "TUBE", idtmed[5], tubpar, 3);
 143   gMC->Gspos("P001", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 144
 145   //-- SECOND SECTION OF THE BEAM PIPE (FROM THE END OF D1 TO THE BEGINNING OF
 146   //    D2)
 147
 148   //-- FROM MAGNETIC BEGINNING OG D1 TO MAGNETIC END OF D1 + 23.5 cm
 149   //-- Elliptic pipe
 150
 151   zd1 = 6310.8-472.5;
 152
 153   elpar[0] = 6.84/2.;
 154   elpar[1] = 5.86/2.;
 155   elpar[2] = 945./2.;
 156   gMC->Gsvolu("E001", "ELTU", idtmed[5], elpar, 3);
 157   gMC->Gspos("E001", 1, "ZDC ", 0., 0., elpar[2] + zd1, 0, "ONLY");
 158
 159   elpar[0] = 6.44/2.;
 160   elpar[1] = 5.46/2.;
 161   elpar[2] = 945./2.;
 162   gMC->Gsvolu("E002", "ELTU", idtmed[10], elpar, 3);
 163   gMC->Gspos("E002", 1, "E001", 0., 0., 0., 0, "ONLY");
 164
 165   zd1 += 2.*elpar[2];
 166
 167   elpar[0] = 6.84/2.;
 168   elpar[1] = 5.86/2.;
 169   elpar[2] = 13.5/2.;
 170   gMC->Gsvolu("E003", "ELTU", idtmed[5], elpar, 3);
 171   gMC->Gspos("E002", 1, "ZDC ", 0., 0., elpar[2] + zd1, 0, "ONLY");
 172
 173   elpar[0] = 6.44/2.;
 174   elpar[1] = 5.46/2.;
 175   elpar[2] = 13.5/2.;
 176   gMC->Gsvolu("E004", "ELTU", idtmed[10], elpar, 3);
 177   gMC->Gspos("E004", 1, "E003", 0., 0., 0., 0, "ONLY");
 178
 179   zd1 += 2.*elpar[2];
 180
 181   conpar[0] = 25./2.;
 182   conpar[1] = 6.44/2.;
 183   conpar[2] = 6.84/2.;
 184   conpar[3] = 10./2.;
 185   conpar[4] = 10.4/2.;
 186   gMC->Gsvolu("C001", "CONE", idtmed[5], conpar, 5);
 187   gMC->Gspos("C001", 1, "ZDC ", 0., 0., conpar[0] + zd1, 0, "ONLY");
 188
 189   zd1 += 2.*conpar[0];
 190
 191   tubpar[0] = 10./2.;
 192   tubpar[1] = 10.4/2.;
 193   tubpar[2] = 50./2.;
 194   gMC->Gsvolu("P002", "TUBE", idtmed[5], tubpar, 3);
 195   gMC->Gspos("P002", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 196
 197   zd1 += tubpar[2] * 2.;
 198
 199   tubpar[0] = 10./2.;
 200   tubpar[1] = 10.4/2.;
 201   tubpar[2] = 10./2.;
 202   gMC->Gsvolu("P003", "TUBE", idtmed[5], tubpar, 3);
 203   gMC->Gspos("P003", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 204
 205   zd1 += tubpar[2] * 2.;
 206
 207   tubpar[0] = 10./2.;
 208   tubpar[1] = 10.4/2.;
 209   tubpar[2] = 3.16/2.;
 210   gMC->Gsvolu("P004", "TUBE", idtmed[5], tubpar, 3);
 211   gMC->Gspos("P004", 1, "ZDC ", 0., 0., tubpar[0] + zd1, 0, "ONLY");
 212
 213   zd1 += tubpar[2] * 2.;
 214
 215   tubpar[0] = 10.0/2.;
 216   tubpar[1] = 10.4/2;
 217   tubpar[2] = 190./2.;
 218   gMC->Gsvolu("P005", "TUBE", idtmed[5], tubpar, 3);
 219   gMC->Gspos("P005", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 220
 221   zd1 += tubpar[2] * 2.;
 222
 223   conpar[0] = 30./2.;
 224   conpar[1] = 10./2.;
 225   conpar[2] = 10.4/2.;
 226   conpar[3] = 20.6/2.;
 227   conpar[4] = 21./2.;
 228   gMC->Gsvolu("P006", "CONE", idtmed[5], conpar, 5);
 229   gMC->Gspos("P006", 1, "ZDC ", 0., 0., conpar[0] + zd1, 0, "ONLY");
 230
 231   zd1 += conpar[0] * 2.;
 232
 233   tubpar[0] = 20.6/2.;
 234   tubpar[1] = 21./2.;
 235   tubpar[2] = 450./2.;
 236   gMC->Gsvolu("P007", "TUBE", idtmed[5], tubpar, 3);
 237   gMC->Gspos("P007", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 238
 239   zd1 += tubpar[2] * 2.;
 240
 241   conpar[0] = 13.6/2.;
 242   conpar[1] = 20.6/2.;
 243   conpar[2] = 21./2.;
 244   conpar[3] = 25.4/2.;
 245   conpar[4] = 25.8/2.;
 246   gMC->Gsvolu("P008", "CONE", idtmed[5], conpar, 5);
 247   gMC->Gspos("P008", 1, "ZDC ", 0., 0., conpar[0] + zd1, 0, "ONLY");
 248
 249   zd1 += conpar[0] * 2.;
 250
 251   tubpar[0] = 25.4/2.;
 252   tubpar[1] = 25.8/2.;
 253   tubpar[2] = 205.8/2.;
 254   gMC->Gsvolu("P009", "TUBE", idtmed[5], tubpar, 3);
 255   gMC->Gspos("P009", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 256
 257   zd1 += tubpar[2] * 2.;
 258
 259   tubpar[0] = 50./2.;
 260   tubpar[1] = 50.4/2.;
 261   tubpar[2] = 505.4/2.;
 262   gMC->Gsvolu("P010", "TUBE", idtmed[5], tubpar, 3);
 263   gMC->Gspos("P010", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 264
 265   zd1 += tubpar[2] * 2.;
 266
 267   tubpar[0] = 50./2.;
 268   tubpar[1] = 50.4/2.;
 269   tubpar[2] = 700./2.;
 270   gMC->Gsvolu("P011", "TUBE", idtmed[5], tubpar, 3);
 271   gMC->Gspos("P011", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 272
 273   zd1 += tubpar[2] * 2.;
 274
 275   tubpar[0] = 50./2.;
 276   tubpar[1] = 50.4/2.;
 277   tubpar[2] = 778.5/2.;
 278   gMC->Gsvolu("P012", "TUBE", idtmed[5], tubpar, 3);
 279   gMC->Gspos("P012", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 280
 281   zd1 += tubpar[2] * 2.;
 282
 283   conpar[0] = 14.18/2.;
 284   conpar[1] = 50./2.;
 285   conpar[2] = 50.4/2.;
 286   conpar[3] = 55./2.;
 287   conpar[4] = 55.4/2.;
 288   gMC->Gsvolu("P013", "CONE", idtmed[5], conpar, 5);
 289   gMC->Gspos("P013", 1, "ZDC ", 0., 0., conpar[0] + zd1, 0, "ONLY");
 290
 291   zd1 += conpar[0] * 2.;
 292
 293   tubpar[0] = 55./2.;
 294   tubpar[1] = 55.4/2.;
 295   tubpar[2] = 730./2.;
 296   gMC->Gsvolu("P014", "TUBE", idtmed[5], tubpar, 3);
 297   gMC->Gspos("P014", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 298
 299   zd1 += tubpar[2] * 2.;
 300
 301   conpar[0] = 36.86/2.;
 302   conpar[1] = 55./2.;
 303   conpar[2] = 55.4/2.;
 304   conpar[3] = 68./2.;
 305   conpar[4] = 68.4/2.;
 306   gMC->Gsvolu("P015", "CONE", idtmed[5], conpar, 5);
 307   gMC->Gspos("P015", 1, "ZDC ", 0., 0., conpar[0] + zd1, 0, "ONLY");
 308
 309   zd1 += conpar[0] * 2.;
 310
 311   tubpar[0] = 68./2.;
 312   tubpar[1] = 68.4/2.;
 313   tubpar[2] = 927.3/2.;
 314   gMC->Gsvolu("P016", "TUBE", idtmed[5], tubpar, 3);
 315   gMC->Gspos("P016", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 316
 317   zd1 += tubpar[2] * 2.;
 318
 319   tubpar[0] = 0./2.;
 320   tubpar[1] = 68.4/2.;
 321   tubpar[2] = 0.2/2.;
 322   gMC->Gsvolu("P017", "TUBE", idtmed[5], tubpar, 3);
 323   gMC->Gspos("P017", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 324
 325   zd1 += tubpar[2] * 2.;
 326
 327   tubpar[0] = 0./2.;
 328   tubpar[1] = 5./2.;
 329   tubpar[2] = 0.2/2.;
 330   gMC->Gsvolu("Q017", "TUBE", idtmed[10], tubpar, 3);
 331
 332   //-- Position Q017 inside P017
 333   gMC->Gspos("Q017", 1, "P017", -7.7, 0., 0., 0, "ONLY");
 334
 335   tubpar[0] = 0./2.;
 336   tubpar[1] = 7./2.;
 337   tubpar[2] = 0.2/2.;
 338   gMC->Gsvolu("R017", "TUBE", idtmed[10], tubpar, 3);
 339
 340   //-- Position R017 inside P017
 341   gMC->Gspos("R017", 1, "P017", 7.7, 0., 0., 0, "ONLY");
 342
 343   //-- BEAM PIPE BETWEEN END OF CONICAL PIPE AND BEGINNING OF D2
 344
 345   tubpar[0] = 5./2.;
 346   tubpar[1] = 5.4/2.;
 347   tubpar[2] = 678./2.;
 348   gMC->Gsvolu("P018", "TUBE", idtmed[5], tubpar, 3);
 349
 350   tubpar[0] = 7./2.;
 351   tubpar[1] = 7.4/2.;
 352   tubpar[2] = 678./2.;
 353   gMC->Gsvolu("P019", "TUBE", idtmed[5], tubpar, 3);
 354
 355   // -- ROTATE PIPES
 356
 357   AliMatrix(im1, 90.-0.071, 0., 90., 90., .071, 180.);
 358   angle = .071*kDegrad;
 359   gMC->Gspos("P018", 1, "ZDC ", TMath::Sin(angle) * 645. / 2. - 9.7 +
 360                TMath::Sin(angle) * 945. / 2., 0., tubpar[2] + zd1, im1, "ONLY");
 361   AliMatrix(im2, 90.+0.071, 0., 90., 90., .071, 0.);
 362   gMC->Gspos("P019", 1, "ZDC ", 9.7 - TMath::Sin(angle) * 645. / 2., 0.,
 363                tubpar[2] + zd1, im2, "ONLY");
 364
 365   // --  END OF BEAM PIPE VOLUME DEFINITION. MAGNET DEFINITION FOLLOWS
 366   //     (LHC OPTICS 6)
 367
 368   // -- COMPENSATOR DIPOLE (MBXW)
 369   //     GAP (VACUUM WITH MAGNETIC FIELD)
 370
 371   tubpar[0] = 0.;
 372   tubpar[1] = 4.5;
 373   tubpar[2] = 340./2.;
 374   gMC->Gsvolu("MBXW", "TUBE", idtmed[11], tubpar, 3);
 375   gMC->Gspos("MBXW", 1, "ZDC ", 0., 0., tubpar[2] + 805., 0, "ONLY");
 376
 377   // --  YOKE (IRON WITHOUT MAGNETIC FIELD)
 378
 379   tubpar[0] = 4.5;
 380   tubpar[1] = 55.;
 381   tubpar[2] = 340./2.;
 382   gMC->Gsvolu("YMBX", "TUBE", idtmed[5], tubpar, 3);
 383   gMC->Gspos("YMBX", 1, "ZDC ", 0., 0., tubpar[2] + 805., 0, "ONLY");
 384
 385   // -- COMPENSATOR DIPOLE (MCBWA)
 386   //     GAP (VACUUM WITH MAGNETIC FIELD)
 387
 388   tubpar[0] = 0.;
 389   tubpar[1] = 4.5;
 390   tubpar[2] = 170./2.;
 391   gMC->Gsvolu("MCBW", "TUBE", idtmed[11], tubpar, 3);
 392   gMC->Gspos("MCBW", 1, "ZDC ", 0., 0., tubpar[2] + 1921.6, 0, "ONLY");
 393
 394   // --  YOKE (IRON WITHOUT MAGNETIC FIELD)
 395
 396   tubpar[0] = 4.5;
 397   tubpar[1] = 55.;
 398   tubpar[2] = 170./2.;
 399   gMC->Gsvolu("YMCB", "TUBE", idtmed[5], tubpar, 3);
 400   gMC->Gspos("YMCB", 1, "ZDC ", 0., 0., tubpar[2] + 1921.6, 0, "ONLY");
 401
 402   // -- INNER TRIPLET
 403
 404   zq = 2300.;
 405
 406   // -- DEFINE MQXL AND MQX QUADRUPOLE ELEMENT
 407
 408   //     MQXL
 409   // --  GAP (VACUUM WITH MAGNETIC FIELD)
 410
 411   tubpar[0] = 0.;
 412   tubpar[1] = 3.5;
 413   tubpar[2] = 630./2.;
 414   gMC->Gsvolu("MQXL", "TUBE", idtmed[11], tubpar, 3);
 415
 416   // --  YOKE
 417
 418   tubpar[0] = 3.5;
 419   tubpar[1] = 22.;
 420   tubpar[2] = 630./2.;
 421   gMC->Gsvolu("YMQL", "TUBE", idtmed[5], tubpar, 3);
 422
 423   gMC->Gspos("MQXL", 1, "ZDC ", 0., 0., tubpar[2] + zq, 0, "ONLY");
 424   gMC->Gspos("YMQL", 1, "ZDC ", 0., 0., tubpar[2] + zq, 0, "ONLY");
 425
 426   gMC->Gspos("MQXL", 2, "ZDC ", 0., 0., tubpar[2] + zq + 2430., 0, "ONLY");
 427   gMC->Gspos("YMQL", 2, "ZDC ", 0., 0., tubpar[2] + zq + 2430., 0, "ONLY");
 428
 429   // --  MQX
 430   // --  GAP (VACUUM WITH MAGNETIC FIELD)
 431
 432   tubpar[0] = 0.;
 433   tubpar[1] = 3.5;
 434   tubpar[2] = 550./2.;
 435   gMC->Gsvolu("MQX ", "TUBE", idtmed[11], tubpar, 3);
 436
 437   // --  YOKE
 438
 439   tubpar[0] = 3.5;
 440   tubpar[1] = 22.;
 441   tubpar[2] = 550./2.;
 442   gMC->Gsvolu("YMQ ", "TUBE", idtmed[5], tubpar, 3);
 443
 444   gMC->Gspos("MQX ", 1, "ZDC ", 0., 0., tubpar[2] + zq + 880.,  0, "ONLY");
 445   gMC->Gspos("YMQ ", 1, "ZDC ", 0., 0., tubpar[2] + zq + 880.,  0, "ONLY");
 446
 447   gMC->Gspos("MQX ", 2, "ZDC ", 0., 0., tubpar[2] + zq + 1530., 0, "ONLY");
 448   gMC->Gspos("YMQ ", 2, "ZDC ", 0., 0., tubpar[2] + zq + 1530., 0, "ONLY");
 449
 450   // -- SEPARATOR DIPOLE D1
 451
 452   zd1 = 5838.3;
 453
 454   // --  GAP (VACUUM WITH MAGNETIC FIELD)
 455
 456   tubpar[0] = 0.;
 457   tubpar[1] = 7.5/2.;
 458   tubpar[2] = 945./2.;
 459   gMC->Gsvolu("D1  ", "TUBE", idtmed[11], tubpar, 3);
 460
 461   // --  YOKE
 462
 463   tubpar[0] = 0.;
 464   tubpar[1] = 110./2;
 465   tubpar[2] = 945./2.;
 466   gMC->Gsvolu("YD1 ", "TUBE", idtmed[5], tubpar, 3);
 467
 468   gMC->Gspos("YD1 ", 1, "ZDC ", 0., 0., tubpar[2] + zd1, 0, "ONLY");
 469   gMC->Gspos("D1  ", 1, "YD1 ", 0., 0., 0., 0, "ONLY");
 470
 471   // -- DIPOLE D2
 472
 473   zd2 = 12147.6;
 474
 475   // --  GAP (VACUUM WITH MAGNETIC FIELD)
 476
 477   tubpar[0] = 0.;
 478   tubpar[1] = 7.5/2.;
 479   tubpar[2] = 945./2.;
 480   gMC->Gsvolu("D2  ", "TUBE", idtmed[11], tubpar, 3);
 481
 482   // --  YOKE
 483
 484   tubpar[0] = 0.;
 485   tubpar[1] = 55.;
 486   tubpar[2] = 945./2.;
 487   gMC->Gsvolu("YD2 ", "TUBE", idtmed[5], tubpar, 3);
 488
 489   gMC->Gspos("YD2 ", 1, "ZDC ", 0., 0., tubpar[2] + zd2, 0, "ONLY");
 490
 491   gMC->Gspos("D2  ", 1, "YD2 ", -9.7, 0., 0., 0, "ONLY");
 492   gMC->Gspos("D2  ", 2, "YD2 ",  9.7, 0., 0., 0, "ONLY");
 493
 494   // -- END OF MAGNET DEFINITION
 495 }
 496
 497 //_____________________________________________________________________________
 498 void AliZDCv1::CreateZDC()
 499 {
 500
 501   Int_t *idtmed = fIdtmed->GetArray();
 502
 503   //-- Create calorimeters geometry
 504
 505   //--> Neutron calorimeter (ZN)
 506
 507   gMC->Gsvolu("ZNEU", "BOX ", idtmed[1], fDimZN, 3); // Passive material
 508   gMC->Gsvolu("ZNF1", "TUBE", idtmed[3], fFibZN, 3); // Active material
 509   gMC->Gsvolu("ZNF2", "TUBE", idtmed[4], fFibZN, 3);
 510   gMC->Gsvolu("ZNF3", "TUBE", idtmed[4], fFibZN, 3);
 511   gMC->Gsvolu("ZNF4", "TUBE", idtmed[3], fFibZN, 3);
 512   gMC->Gsvolu("ZNG1", "BOX ", idtmed[12], fGrvZN, 3); // Empty grooves
 513   gMC->Gsvolu("ZNG2", "BOX ", idtmed[12], fGrvZN, 3);
 514   gMC->Gsvolu("ZNG3", "BOX ", idtmed[12], fGrvZN, 3);
 515   gMC->Gsvolu("ZNG4", "BOX ", idtmed[12], fGrvZN, 3);
 516
 517   // Divide ZNEU in towers (for hits purposes)
 518
 519   gMC->Gsdvn("ZNTX", "ZNEU", fTowZN[0], 1); // x-tower
 520   gMC->Gsdvn("ZN1 ", "ZNTX", fTowZN[1], 2); // y-tower
 521
 522   //-- Divide ZN1 in minitowers
 523   //  fDivZN[0]= NUMBER OF FIBERS PER TOWER ALONG X-AXIS,
 524   //  fDivZN[1]= NUMBER OF FIBERS PER TOWER ALONG Y-AXIS
 525   //  (4 fibres per minitower)
 526
 527   gMC->Gsdvn("ZNSL", "ZN1 ", fDivZN[1], 2); // Slices
 528   gMC->Gsdvn("ZNST", "ZNSL", fDivZN[0], 1); // Sticks
 529
 530   // --- Position the empty grooves in the sticks (4 grooves per stick)
 531   Float_t dx = fDimZN[0] / fDivZN[0] / 4.;
 532   Float_t dy = fDimZN[1] / fDivZN[1] / 4.;
 533
 534   gMC->Gspos("ZNG1", 1, "ZNST", 0.-dx, 0.+dy, 0., 0, "ONLY");
 535   gMC->Gspos("ZNG2", 1, "ZNST", 0.+dx, 0.+dy, 0., 0, "ONLY");
 536   gMC->Gspos("ZNG3", 1, "ZNST", 0.-dx, 0.-dy, 0., 0, "ONLY");
 537   gMC->Gspos("ZNG4", 1, "ZNST", 0.+dx, 0.-dy, 0., 0, "ONLY");
 538
 539   // --- Position the fibers in the grooves
 540   gMC->Gspos("ZNF1", 1, "ZNG1", 0., 0., 0., 0, "ONLY");
 541   gMC->Gspos("ZNF2", 1, "ZNG2", 0., 0., 0., 0, "ONLY");
 542   gMC->Gspos("ZNF3", 1, "ZNG3", 0., 0., 0., 0, "ONLY");
 543   gMC->Gspos("ZNF4", 1, "ZNG4", 0., 0., 0., 0, "ONLY");
 544
 545   // --- Position the neutron calorimeter in ZDC
 546   gMC->Gspos("ZNEU", 1, "ZDC ", fPosZN[0], fPosZN[1], fPosZN[2] + fDimZN[2], 0, "ONLY");
 547
 548   //--> Proton calorimeter
 549
 550   gMC->Gsvolu("ZPRO", "BOX ", idtmed[2], fDimZP, 3); // Passive material
 551   gMC->Gsvolu("ZPF1", "TUBE", idtmed[3], fFibZP, 3); // Active material
 552   gMC->Gsvolu("ZPF2", "TUBE", idtmed[4], fFibZP, 3);
 553   gMC->Gsvolu("ZPF3", "TUBE", idtmed[4], fFibZP, 3);
 554   gMC->Gsvolu("ZPF4", "TUBE", idtmed[3], fFibZP, 3);
 555   gMC->Gsvolu("ZPG1", "BOX ", idtmed[12], fGrvZP, 3); // Empty grooves
 556   gMC->Gsvolu("ZPG2", "BOX ", idtmed[12], fGrvZP, 3);
 557   gMC->Gsvolu("ZPG3", "BOX ", idtmed[12], fGrvZP, 3);
 558   gMC->Gsvolu("ZPG4", "BOX ", idtmed[12], fGrvZP, 3);
 559
 560   //-- Divide ZPRO in towers(for hits purposes)
 561
 562   gMC->Gsdvn("ZPTX", "ZPRO", fTowZP[0], 1); // x-tower
 563   gMC->Gsdvn("ZP1 ", "ZPTX", fTowZP[1], 2); // y-tower
 564
 565
 566   //-- Divide ZP1 in minitowers
 567   //  fDivZP[0]= NUMBER OF FIBERS ALONG X-AXIS PER MINITOWER,
 568   //  fDivZP[1]= NUMBER OF FIBERS ALONG Y-AXIS PER MINITOWER
 569   //  (4 fiber per minitower)
 570
 571   gMC->Gsdvn("ZPSL", "ZP1 ", fDivZP[1], 2); // Slices
 572   gMC->Gsdvn("ZPST", "ZPSL", fDivZP[0], 1); // Sticks
 573
 574   // --- Position the empty grooves in the sticks (4 grooves per stick)
 575   dx = fDimZP[0] / fTowZP[0] / fDivZP[0] / 2.;
 576   dy = fDimZP[1] / fTowZP[1] / fDivZP[1] / 2.;
 577
 578   gMC->Gspos("ZPG1", 1, "ZPST", 0.-dx, 0.+dy, 0., 0, "ONLY");
 579   gMC->Gspos("ZPG2", 1, "ZPST", 0.+dx, 0.+dy, 0., 0, "ONLY");
 580   gMC->Gspos("ZPG3", 1, "ZPST", 0.-dx, 0.-dy, 0., 0, "ONLY");
 581   gMC->Gspos("ZPG4", 1, "ZPST", 0.+dx, 0.-dy, 0., 0, "ONLY");
 582
 583   // --- Position the fibers in the grooves
 584   gMC->Gspos("ZPF1", 1, "ZPG1", 0., 0., 0., 0, "ONLY");
 585   gMC->Gspos("ZPF2", 1, "ZPG2", 0., 0., 0., 0, "ONLY");
 586   gMC->Gspos("ZPF3", 1, "ZPG3", 0., 0., 0., 0, "ONLY");
 587   gMC->Gspos("ZPF4", 1, "ZPG4", 0., 0., 0., 0, "ONLY");
 588
 589
 590   // --- Position the proton calorimeter in ZDC
 591   gMC->Gspos("ZPRO", 1, "ZDC ", fPosZP[0], fPosZP[1], fPosZP[2] + fDimZP[2], 0, "ONLY");
 592
 593 }
 594
 595 //_____________________________________________________________________________
 596 void AliZDCv1::DrawModule()
 597 {
 598   //
 599   // Draw a shaded view of the Zero Degree Calorimeter version 1
 600   //
 601
 602   // Set everything unseen
 603   gMC->Gsatt("*", "seen", -1);
 604   //
 605   // Set ALIC mother transparent
 606   gMC->Gsatt("ALIC","SEEN",0);
 607   //
 608   // Set the volumes visible
 609   gMC->Gsatt("ZDC ","SEEN",0);
 610   gMC->Gsatt("P001","SEEN",1);
 611   gMC->Gsatt("E001","SEEN",1);
 612   gMC->Gsatt("E002","SEEN",1);
 613   gMC->Gsatt("E003","SEEN",1);
 614   gMC->Gsatt("E004","SEEN",1);
 615   gMC->Gsatt("C001","SEEN",1);
 616   gMC->Gsatt("P002","SEEN",1);
 617   gMC->Gsatt("P003","SEEN",1);
 618   gMC->Gsatt("P004","SEEN",1);
 619   gMC->Gsatt("P005","SEEN",1);
 620   gMC->Gsatt("P006","SEEN",1);
 621   gMC->Gsatt("P007","SEEN",1);
 622   gMC->Gsatt("P008","SEEN",1);
 623   gMC->Gsatt("P009","SEEN",1);
 624   gMC->Gsatt("P010","SEEN",1);
 625   gMC->Gsatt("P011","SEEN",1);
 626   gMC->Gsatt("P012","SEEN",1);
 627   gMC->Gsatt("P013","SEEN",1);
 628   gMC->Gsatt("P014","SEEN",1);
 629   gMC->Gsatt("P015","SEEN",1);
 630   gMC->Gsatt("P016","SEEN",1);
 631   gMC->Gsatt("P017","SEEN",1);
 632   gMC->Gsatt("Q017","SEEN",1);
 633   gMC->Gsatt("R017","SEEN",1);
 634   gMC->Gsatt("P018","SEEN",1);
 635   gMC->Gsatt("P019","SEEN",1);
 636   gMC->Gsatt("MBXW","SEEN",1);
 637   gMC->Gsatt("YMBX","SEEN",1);
 638   gMC->Gsatt("MCBW","SEEN",1);
 639   gMC->Gsatt("YMCB","SEEN",1);
 640   gMC->Gsatt("MQXL","SEEN",1);
 641   gMC->Gsatt("YMQL","SEEN",1);
 642   gMC->Gsatt("MQX ","SEEN",1);
 643   gMC->Gsatt("YMQ ","SEEN",1);
 644   gMC->Gsatt("D1  ","SEEN",1);
 645   gMC->Gsatt("YD1 ","SEEN",1);
 646   gMC->Gsatt("D2  ","SEEN",1);
 647   gMC->Gsatt("YD2 ","SEEN",1);
 648   gMC->Gsatt("ZNEU","SEEN",0);
 649   gMC->Gsatt("ZNF1","SEEN",0);
 650   gMC->Gsatt("ZNF2","SEEN",0);
 651   gMC->Gsatt("ZNF3","SEEN",0);
 652   gMC->Gsatt("ZNF4","SEEN",0);
 653   gMC->Gsatt("ZNG1","SEEN",0);
 654   gMC->Gsatt("ZNG2","SEEN",0);
 655   gMC->Gsatt("ZNG3","SEEN",0);
 656   gMC->Gsatt("ZNG4","SEEN",0);
 657   gMC->Gsatt("ZNTX","SEEN",0);
 658   gMC->Gsatt("ZN1 ","COLO",2);
 659   gMC->Gsatt("ZN1 ","SEEN",1);
 660   gMC->Gsatt("ZNSL","SEEN",0);
 661   gMC->Gsatt("ZNST","SEEN",0);
 662   gMC->Gsatt("ZPRO","SEEN",0);
 663   gMC->Gsatt("ZPF1","SEEN",0);
 664   gMC->Gsatt("ZPF2","SEEN",0);
 665   gMC->Gsatt("ZPF3","SEEN",0);
 666   gMC->Gsatt("ZPF4","SEEN",0);
 667   gMC->Gsatt("ZPG1","SEEN",0);
 668   gMC->Gsatt("ZPG2","SEEN",0);
 669   gMC->Gsatt("ZPG3","SEEN",0);
 670   gMC->Gsatt("ZPG4","SEEN",0);
 671   gMC->Gsatt("ZPTX","SEEN",0);
 672   gMC->Gsatt("ZP1 ","COLO",2);
 673   gMC->Gsatt("ZP1 ","SEEN",1);
 674   gMC->Gsatt("ZPSL","SEEN",0);
 675   gMC->Gsatt("ZPST","SEEN",0);
 676
 677   //
 678   gMC->Gdopt("hide", "on");
 679   gMC->Gdopt("shad", "on");
 680   gMC->Gsatt("*", "fill", 7);
 681   gMC->SetClipBox(".");
 682   gMC->SetClipBox("*", 0, 100, -100, 100, 12000, 16000);
 683   gMC->DefaultRange();
 684   gMC->Gdraw("alic", 40, 30, 0, 488, 220, .07, .07);
 685   gMC->Gdhead(1111, "Zero Degree Calorimeter Version 1");
 686   gMC->Gdman(18, 4, "MAN");
 687 }
 688
 689 //_____________________________________________________________________________
 690 void AliZDCv1::CreateMaterials()
 691 {
 692   //
 693   // Create Materials for the Zero Degree Calorimeter
 694   //
 695   // Origin    : E. Scomparin
 696
 697   Int_t *idtmed = fIdtmed->GetArray();
 698
 699   Float_t dens, ubuf[1], wmat[2];
 700   Int_t isvolActive;
 701   Float_t a[2];
 702   Int_t i;
 703   Float_t z[2], epsil=0.001, stmin=0.01;
 704   Int_t isvol;
 705   Float_t fieldm = gAlice->Field()->Max();
 706   Int_t inofld;
 707   Float_t deemax=-1;
 708   Float_t tmaxfd=gAlice->Field()->Max();
 709   Int_t isxfld = gAlice->Field()->Integ();
 710   Float_t stemax;
 711
 712   // --- Store in UBUF r0 for nuclear radius calculation R=r0*A**1/3
 713
 714   // --- Tantalum -> ZN passive material
 715   ubuf[0] = 1.1;
 716   AliMaterial(1, "TANT", 180.95, 73., 16.65, .4, 11.9, ubuf, 1);
 717
 718   // --- Tungsten
 719 //  ubuf[0] = 1.11;
 720 //  AliMaterial(1, "TUNG", 183.85, 74., 19.3, .35, 10.3, ubuf, 1);
 721
 722   // --- Brass (CuZn)  -> ZP passive material
 723   dens = 8.48;
 724   a[0] = 63.546;
 725   a[1] = 65.39;
 726   z[0] = 29.;
 727   z[1] = 30.;
 728   wmat[0] = .63;
 729   wmat[1] = .37;
 730   AliMixture(2, "BRASS               ", a, z, dens, 2, wmat);
 731
 732   // --- SiO2
 733   dens = 2.64;
 734   a[0] = 28.086;
 735   a[1] = 15.9994;
 736   z[0] = 14.;
 737   z[1] = 8.;
 738   wmat[0] = 1.;
 739   wmat[1] = 2.;
 740   AliMixture(3, "SIO2                ", a, z, dens, -2, wmat);
 741
 742   // --- Copper
 743 //  ubuf[0] = 1.1;
 744 //  AliMaterial(7, "COPP", 63.54, 29., 8.96, 1.4, 0., ubuf, 1);
 745
 746
 747   // --- Lead
 748 //  ubuf[0] = 1.12;
 749 //  AliMaterial(6, "LEAD", 207.19, 82., 11.35, .56, 18.5, ubuf, 1);
 750
 751   // --- Iron
 752   ubuf[0] = 1.1;
 753   AliMaterial(5, "IRON", 55.85, 26., 7.87, 1.76, 0., ubuf, 1);
 754
 755   // --- Vacuum (no magnetic field)
 756   AliMaterial(10, "VOID", 1e-16, 1e-16, 1e-16, 1e16, 1e16, ubuf,0);
 757
 758   // --- Vacuum (with magnetic field)
 759   AliMaterial(11, "VOIM", 1e-16, 1e-16, 1e-16, 1e16, 1e16, ubuf,0);
 760
 761   // --- Air (no magnetic field)
 762   AliMaterial(12, "Air    $", 14.61, 7.3, .001205, 30420., 67500., ubuf, 0);
 763
 764   // ---  Definition of tracking media:
 765
 766   // --- Tantalum = 1 ;
 767   // --- Brass = 2 ;
 768   // --- Fibers (SiO2) = 3 ;
 769   // --- Fibers (SiO2) = 4 ;
 770   // --- Iron = 5 ;
 771   // --- Lead = 6 ;
 772   // --- Vacuum (no field) = 10
 773   // --- Vacuum (with field) = 11
 774   // --- Air (no field) = 12
 775
 776
 777   // --- Tracking media parameters
 778   epsil  = .01;
 779   stemax = 1.;
 780   isvol  = 0;
 781   isvolActive = 1;
 782   inofld = 0;
 783   fieldm = 0.;
 784
 785   AliMedium(1, "ZTANT", 1, isvolActive, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 786 //  AliMedium(1, "ZW", 1, isvolActive, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 787   AliMedium(2, "ZBRASS", 2, isvolActive, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 788   AliMedium(3, "ZSIO2", 3, isvolActive, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 789   AliMedium(4, "ZQUAR", 3, isvolActive, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 790 //  AliMedium(7, "ZCOPP", 7, isvolActive, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 791 //  AliMedium(6, "ZLEAD", 6, isvolActive, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 792   AliMedium(5, "ZIRON", 5, isvol, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 793   AliMedium(10, "ZVOID", 10, isvol, inofld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 794   AliMedium(12, "ZAIR", 12, 0, inofld, fieldm, tmaxfd, stemax,deemax, epsil, stmin);
 795
 796   fieldm = 45.;
 797   AliMedium(11, "ZVOIM", 11, isvol, isxfld, fieldm, tmaxfd, stemax, deemax, epsil, stmin);
 798
 799   // Thresholds for showering in the ZDCs
 800
 801   i = 1;
 802   gMC->Gstpar(idtmed[i], "CUTGAM", .001);
 803   gMC->Gstpar(idtmed[i], "CUTELE", .001);
 804   gMC->Gstpar(idtmed[i], "CUTNEU", .01);
 805   gMC->Gstpar(idtmed[i], "CUTHAD", .01);
 806   i = 2;
 807   gMC->Gstpar(idtmed[i], "CUTGAM", .001);
 808   gMC->Gstpar(idtmed[i], "CUTELE", .001);
 809   gMC->Gstpar(idtmed[i], "CUTNEU", .01);
 810   gMC->Gstpar(idtmed[i], "CUTHAD", .01);
 811
 812   // Avoid too detailed showering along the beam line
 813
 814   i = 5;
 815   gMC->Gstpar(idtmed[i], "CUTGAM", .1);
 816   gMC->Gstpar(idtmed[i], "CUTELE", .1);
 817   gMC->Gstpar(idtmed[i], "CUTNEU", 1.);
 818   gMC->Gstpar(idtmed[i], "CUTHAD", 1.);
 819
 820   // Avoid interaction in fibers (only energy loss allowed)
 821   i = 3;
 822   gMC->Gstpar(idtmed[i], "DCAY", 0.);
 823   gMC->Gstpar(idtmed[i], "MULS", 0.);
 824   gMC->Gstpar(idtmed[i], "PFIS", 0.);
 825   gMC->Gstpar(idtmed[i], "MUNU", 0.);
 826   gMC->Gstpar(idtmed[i], "LOSS", 1.);
 827   gMC->Gstpar(idtmed[i], "PHOT", 0.);
 828   gMC->Gstpar(idtmed[i], "COMP", 0.);
 829   gMC->Gstpar(idtmed[i], "PAIR", 0.);
 830   gMC->Gstpar(idtmed[i], "BREM", 0.);
 831   gMC->Gstpar(idtmed[i], "DRAY", 0.);
 832   gMC->Gstpar(idtmed[i], "ANNI", 0.);
 833   gMC->Gstpar(idtmed[i], "HADR", 0.);
 834   i = 4;
 835   gMC->Gstpar(idtmed[i], "DCAY", 0.);
 836   gMC->Gstpar(idtmed[i], "MULS", 0.);
 837   gMC->Gstpar(idtmed[i], "PFIS", 0.);
 838   gMC->Gstpar(idtmed[i], "MUNU", 0.);
 839   gMC->Gstpar(idtmed[i], "LOSS", 1.);
 840   gMC->Gstpar(idtmed[i], "PHOT", 0.);
 841   gMC->Gstpar(idtmed[i], "COMP", 0.);
 842   gMC->Gstpar(idtmed[i], "PAIR", 0.);
 843   gMC->Gstpar(idtmed[i], "BREM", 0.);
 844   gMC->Gstpar(idtmed[i], "DRAY", 0.);
 845   gMC->Gstpar(idtmed[i], "ANNI", 0.);
 846   gMC->Gstpar(idtmed[i], "HADR", 0.);
 847   //
 848   fMedSensF1 = idtmed[3];  // Sensitive volume: fibres type 1
 849   fMedSensF2 = idtmed[4];  // Sensitive volume: fibres type 2
 850   fMedSensZN = idtmed[1];  // Sensitive volume: ZN passive material
 851   fMedSensZP = idtmed[2];  // Sensitive volume: ZP passive material
 852   fMedSensGR = idtmed[12]; // Sensitive volume: air into the grooves
 853 }
 854
 855 //_____________________________________________________________________________
 856 void AliZDCv1::Init()
 857 {
 858  InitTables();
 859
 860 }
 861
 862 //_____________________________________________________________________________
 863 void AliZDCv1::InitTables()
 864 {
 865   Int_t k, j;
 866   //Initialize parameters for light tables and read them
 867   fNalfan = 90;
 868   fNalfap = 90;
 869   fNben = 18;
 870   fNbep = 28;
 871
 872   FILE *fp1, *fp2, *fp3, *fp4, *fp5, *fp6, *fp7, *fp8;
 873
 874   if((fp1 = fopen("light22620362207s","r")) == NULL){
 875      printf("Cannot open file fp1 \n");
 876      return;
 877   }
 878   if((fp2 = fopen("light22620362208s","r")) == NULL){
 879      printf("Cannot open file fp2 \n");
 880      return;
 881   }
 882   if((fp3 = fopen("light22620362209s","r")) == NULL){
 883      printf("Cannot open file fp3 \n");
 884      return;
 885   }
 886   if((fp4 = fopen("light22620362210s","r")) == NULL){
 887      printf("Cannot open file fp4 \n");
 888      return;
 889   }
 890 //  printf(" --- Reading light tables for ZN \n");
 891   for(k=0; k<fNalfan; k++){
 892      for(j=0; j<fNben; j++){
 893        fscanf(fp1,"%f",&fTablen[0][k][j]);
 894        fscanf(fp2,"%f",&fTablen[1][k][j]);
 895        fscanf(fp3,"%f",&fTablen[2][k][j]);
 896        fscanf(fp4,"%f",&fTablen[3][k][j]);
 897      }
 898   }
 899   fclose(fp1);
 900   fclose(fp2);
 901   fclose(fp3);
 902   fclose(fp4);
 903
 904   if((fp5 = fopen("light22620552207s","r")) == NULL){
 905      printf("Cannot open file fp5 \n");
 906      return;
 907   }
 908   if((fp6 = fopen("light22620552208s","r")) == NULL){
 909      printf("Cannot open file fp6 \n");
 910      return;
 911   }
 912   if((fp7 = fopen("light22620552209s","r")) == NULL){
 913      printf("Cannot open file fp7 \n");
 914      return;
 915   }
 916   if((fp8 = fopen("light22620552210s","r")) == NULL){
 917      printf("Cannot open file fp8 \n");
 918      return;
 919   }
 920 //  printf(" --- Reading light tables for ZP \n");
 921   for(k=0; k<fNalfap; k++){
 922      for(j=0; j<fNbep; j++){
 923        fscanf(fp5,"%f",&fTablep[0][k][j]);
 924        fscanf(fp6,"%f",&fTablep[1][k][j]);
 925        fscanf(fp7,"%f",&fTablep[2][k][j]);
 926        fscanf(fp8,"%f",&fTablep[3][k][j]);
 927      }
 928   }
 929   fclose(fp5);
 930   fclose(fp6);
 931   fclose(fp7);
 932   fclose(fp8);
 933 }
 934
 935 //_____________________________________________________________________________
 936 void AliZDCv1::StepManager()
 937 {
 938   //
 939   // Routine called at every step in the Zero Degree Calorimeters
 940   //
 941
 942   Int_t j;
 943
 944   Int_t vol[2], ibeta, ialfa, ibe;
 945   Float_t x[3], xdet[3], destep, hits[10], m, ekin, um[3], ud[3], be, radius, out;
 946   TLorentzVector s, p;
 947   const char *knamed;
 948
 949   if((gMC->GetMedium() == fMedSensZN) || (gMC->GetMedium() == fMedSensZP) ||
 950      (gMC->GetMedium() == fMedSensGR) || (gMC->GetMedium() == fMedSensF1) ||
 951      (gMC->GetMedium() == fMedSensF2)){
 952
 953   //Particle coordinates
 954     gMC->TrackPosition(s);
 955     for(j=0; j<=2; j++){
 956        x[j] = s[j];
 957     }
 958     hits[0] = x[0];
 959     hits[1] = x[1];
 960     hits[2] = x[2];
 961
 962   // Determine in which ZDC the particle is
 963     knamed = gMC->CurrentVolName();
 964     if(!strncmp(knamed,"ZN",2))vol[0]=1;
 965     if(!strncmp(knamed,"ZP",2))vol[0]=2;
 966
 967   // Determine in which quadrant the particle is
 968
 969     //Quadrant in ZN
 970     gMC->Gmtod(x,xdet,1);
 971     if(vol[0]==1){
 972       if((xdet[0]<0.) && (xdet[1]>0.)) vol[1]=1;
 973       if((xdet[0]>0.) && (xdet[1]>0.)) vol[1]=2;
 974       if((xdet[0]<0.) && (xdet[1]<0.)) vol[1]=3;
 975       if((xdet[0]>0.) && (xdet[1]<0.)) vol[1]=4;
 976     }
 977
 978     //Quadrant in ZP
 979     if(vol[0]==2){
 980       Float_t xqZP = xdet[0]/(fDimZP[0]/2);
 981       for(int i=1; i<=4; i++){
 982          if(xqZP>(i-3) && xqZP<(i-2)){
 983            vol[1] = i;
 984            break;
 985         }
 986      }
 987     }
 988 //    printf("  -> Det. %d Quad. %d \n", vol[0], vol[1]);
 989
 990   // Store impact point and kinetic energy of the ENTERING particle
 991
 992 //    Int_t Curtrack = gAlice->CurrentTrack();
 993 //    Int_t Prim = gAlice->GetPrimary(Curtrack);
 994 //    printf ("Primary: %d, Current Track: %d \n", Prim, Curtrack);
 995
 996 //    if(Curtrack==Prim){
 997       if(gMC->IsTrackEntering()){
 998         //Particle energy
 999         gMC->TrackMomentum(p);
1000 //       printf("p[0] = %f, p[1] = %f, p[2] = %f, p[3] = %f \n",
1001 //                 p[0], p[1], p[2], p[3]);
1002         hits[3] = p[3];
1003
1004         // Impact point on ZN
1005         hits[4] = xdet[0];
1006         hits[5] = xdet[1];
1007         hits[7] = 0;
1008         hits[8] = 0;
1009         hits[9] = 0;
1010
1011 //        printf(" hits[2] = %f \n",hits[2]);
1012         AddHit(gAlice->CurrentTrack(), vol, hits);
1013       }
1014 //    }
1015
1016       // Charged particles -> Energy loss
1017       if((destep=gMC->Edep())){
1018          if(gMC->IsTrackStop()){
1019            gMC->TrackMomentum(p);
1020            m = gMC->TrackMass();
1021            ekin = p[3]-m;
1022            if(ekin<0.) printf("ATTENTION!!!!!!!!!!!!!!! ->      ekin = %f <0 (?)",ekin);
1023            hits[9] = ekin;
1024            hits[7] = 0.;
1025            hits[8] = 0.;
1026            AddHit(gAlice->CurrentTrack(), vol, hits);
1027            }
1028          else{
1029            hits[9] = destep;
1030            hits[7] = 0.;
1031            hits[8] = 0.;
1032            AddHit(gAlice->CurrentTrack(), vol, hits);
1033            }
1034 //       printf("       -> Charged particle -> Dep. E = %f eV \n",hits[8]);
1035          }
1036 //       printf(" \n");
1037   }
1038
1039
1040   // *** Light production in fibres
1041   if((gMC->GetMedium() == fMedSensF1) || (gMC->GetMedium() == fMedSensF2)){
1042 //    printf("%%%%%%%%%%%%%%%% Particle in fibre %%%%%%%%%%%%%%%%\n");
1043
1044      //Select charged particles
1045      if((destep=gMC->Edep())){
1046 //       printf("               -> CHARGED particle!!! \n");
1047
1048        // Particle velocity
1049        gMC->TrackMomentum(p);
1050        Float_t ptot=TMath::Sqrt(p[0]*p[0]+p[1]*p[1]+p[2]*p[2]);
1051        Float_t beta =  ptot/p[3];
1052 //       printf("p[0] = %f, p[1] = %f, p[2] = %f, p[3] = %f, ptot = %f \n",
1053 //              p[0], p[1], p[2], p[3], ptot);
1054 //       Int_t pcID = gMC->TrackPid();
1055 //       printf("       Pc %d in quadrant %d -> beta = %f \n", pcID, vol[1], beta);
1056        if(beta<0.67) return;
1057        else if((beta>=0.67) && (beta<=0.75)) ibeta = 0;
1058        else if((beta>0.75)  && (beta<=0.85)) ibeta = 1;
1059        else if((beta>0.85)  && (beta<=0.95)) ibeta = 2;
1060        else ibeta = 3;
1061
1062        // Angle between particle trajectory and fibre axis
1063        // 1 -> Momentum directions
1064        um[0] = p[0]/ptot;
1065        um[1] = p[1]/ptot;
1066        um[2] = p[2]/ptot;
1067        gMC->Gmtod(um,ud,2);
1068        // 2 -> Angle < limit angle
1069        Double_t alfar = TMath::ACos(ud[2]);
1070        Double_t alfa = alfar*kRaddeg;
1071        if(alfa>110.) return;
1072        ialfa = Int_t(1.+alfa/2.);
1073
1074        // Distance between particle trajectory and fibre axis
1075        gMC->TrackPosition(s);
1076        for(j=0; j<=2; j++){
1077           x[j] = s[j];
1078        }
1079        gMC->Gmtod(x,xdet,1);
1080        if(TMath::Abs(ud[0])>0.00001){
1081          Float_t dcoeff = ud[1]/ud[0];
1082          be = TMath::Abs((xdet[1]-dcoeff*xdet[0])/TMath::Sqrt(dcoeff*dcoeff+1.));
1083        }
1084        else{
1085          be = TMath::Abs(ud[0]);
1086        }
1087
1088        if((vol[0]==1)) radius = fFibZN[1];
1089        if((vol[0]==2)) radius = fFibZP[1];
1090        ibe = Int_t(be*1000.+1);
1091
1092        //Looking into the light tables
1093        Float_t charge = gMC->TrackCharge();
1094
1095        // (1)  ZN
1096        if((vol[0]==1)) {
1097          if(ibe>fNben) ibe=fNben;
1098          out =  charge*charge*fTablen[ibeta][ialfa][ibe];
1099 //       printf("       -> fTablen [%d][%d][%d] = %f \n",
1100 //              ibeta, ialfa, ibe, fTablen[ibeta][ialfa][ibe]);
1101          if(gMC->GetMedium() == fMedSensF1){
1102            hits[7] = out;       //fLightPMQ
1103            hits[8] = 0;
1104            hits[9] = 0;
1105            AddHit(gAlice->CurrentTrack(), vol, hits);
1106          }
1107          else{
1108            hits[7] = 0;
1109            hits[8] = out;       //fLightPMC
1110            hits[9] = 0;
1111            AddHit(gAlice->CurrentTrack(), vol, hits);
1112          }
1113        }
1114
1115        // (2) ZP
1116        if((vol[0]==2)) {
1117          if(ibe>fNbep) ibe=fNbep;
1118          out =  charge*charge*fTablep[ibeta][ialfa][ibe];
1119 //       printf("       -> fTablep [%d][%d][%d] = %f\n",
1120 //              ibeta, ialfa, ibe, fTablen[ibeta][ialfa][ibe]);
1121          if(gMC->GetMedium() == fMedSensF1){
1122            hits[7] = out;       //fLightPMQ
1123            hits[8] = 0;
1124            hits[9] = 0;
1125            AddHit(gAlice->CurrentTrack(), vol, hits);
1126          }
1127          else{
1128            hits[7] = 0;
1129            hits[8] = out;       //fLightPMC
1130            hits[9] = 0;
1131            AddHit(gAlice->CurrentTrack(), vol, hits);
1132          }
1133        }
1134      }
1135 //    printf("\n");
1136    }
1137 }