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