[u/mrichter/AliRoot.git] / FMD / AliFMDSimulator.cxx

/**************************************************************************
 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
 *                                                                        *
 * Author: The ALICE Off-line Project.                                    *
 * Contributors are mentioned in the code where appropriate.              *
 *                                                                        *
 * Permission to use, copy, modify and distribute this software and its   *
 * documentation strictly for non-commercial purposes is hereby granted   *
 * without fee, provided that the above copyright notice appears in all   *
 * copies and that both the copyright notice and this permission notice   *
 * appear in the supporting documentation. The authors make no claims     *
 * about the suitability of this software for any purpose. It is          *
 * provided "as is" without express or implied warranty.                  *
 **************************************************************************/

/* $Id$ */

//____________________________________________________________________
//                                                                          
// Forward Multiplicity Detector based on Silicon wafers. This class
// contains the base procedures for the Forward Multiplicity detector
// Detector consists of 3 sub-detectors FMD1, FMD2, and FMD3, each of
// which has 1 or 2 rings of silicon sensors. 
//                                                       
// This is the base class for all FMD manager classes. 
//                    
// The actual code is done by various separate classes.   Below is
// diagram showing the relationship between the various FMD classes
// that handles the simulation
//
//      +--------+ 1     +-----------------+ 
//      | AliFMD |<>-----| AliFMDSimulator |
//      +--------+	 +-----------------+
//                               ^              
//                               |
//                 +-------------+-------------+
//                 |                           |	      
//        +--------------------+   +-------------------+
//        | AliFMDGeoSimulator |   | AliFMDG3Simulator | 
//        +--------------------+   +---------+---------+
//      
//
// *  AliFMD 
//    This defines the interface for the various parts of AliROOT that
//    uses the FMD, like AliFMDSimulator, AliFMDDigitizer, 
//    AliFMDReconstructor, and so on. 
//
// *  AliFMDSimulator
//    This is the base class for the FMD simulation tasks.   The
//    simulator tasks are responsible to implment the geoemtry, and
//    process hits. 
//                                                                          
// *  AliFMDGeoSimulator
//    This is a concrete implementation of the AliFMDSimulator that
//    uses the TGeo classes directly only.  This defines the active
//    volume as an ONLY XTRU shape with a divided MANY TUBS shape
//    inside to implement the particular shape of the silicon
//    sensors. 
//
// *  AliFMDG3Simulator
//    This is a concrete implementation of the AliFMDSimulator that
//    uses the TVirtualMC interface with GEANT 3.21-like messages.
//    This implements the active volume as a divided TUBS shape.  Hits
//    in the corners should be cut away at run time (but currently
//    isn't). 
//
#include "AliFMDSimulator.h"	// ALIFMDSIMULATOR_H
#include "AliFMDGeometry.h"	// ALIFMDGEOMETRY_H
#include "AliFMDDetector.h"	// ALIFMDDETECTOR_H
#include "AliFMDRing.h"		// ALIFMDRING_H
#include "AliFMD1.h"		// ALIFMD1_H
#include "AliFMD2.h"		// ALIFMD2_H
#include "AliFMD3.h"		// ALIFMD3_H
#include "AliFMD.h"		// ALIFMD_H
#include <AliRun.h>		// ALIRUN_H
#include <AliMC.h>		// ALIMC_H
#include <AliMagF.h>		// ALIMAGF_H
#include <AliLog.h>		// ALILOG_H
#include <TGeoVolume.h>		// ROOT_TGeoVolume
#include <TGeoTube.h>		// ROOT_TGeoTube
#include <TGeoPcon.h>		// ROOT_TGeoPcon
#include <TGeoMaterial.h>	// ROOT_TGeoMaterial
#include <TGeoMedium.h>		// ROOT_TGeoMedium
#include <TGeoXtru.h>		// ROOT_TGeoXtru
#include <TGeoPolygon.h>	// ROOT_TGeoPolygon
#include <TGeoTube.h>		// ROOT_TGeoTube
#include <TGeoManager.h>	// ROOT_TGeoManager
#include <TTree.h>		// ROOT_TTree
#include <TParticle.h>		// ROOT_TParticle
#include <TLorentzVector.h>	// ROOT_TLorentzVector
#include <TVector2.h>		// ROOT_TVector2
#include <TVector3.h>		// ROOT_TVector3
#include <TVirtualMC.h>		// ROOT_TVirtualMC
#include <TArrayD.h>		// ROOT_TArrayD

//====================================================================
ClassImp(AliFMDSimulator)
#if 0
  ; // This is here to keep Emacs for indenting the next line
#endif

//____________________________________________________________________
const Char_t* AliFMDSimulator::fgkActiveName	= "F%cAC";
const Char_t* AliFMDSimulator::fgkSectorName	= "F%cSE";
const Char_t* AliFMDSimulator::fgkStripName	= "F%cST";
const Char_t* AliFMDSimulator::fgkModuleName	= "F%cMO";
const Char_t* AliFMDSimulator::fgkPCBName	= "F%cP%c";
const Char_t* AliFMDSimulator::fgkLongLegName	= "F%cLL";
const Char_t* AliFMDSimulator::fgkShortLegName	= "F%cSL";
const Char_t* AliFMDSimulator::fgkFrontVName	= "F%cFV";
const Char_t* AliFMDSimulator::fgkBackVName	= "F%cBV";
const Char_t* AliFMDSimulator::fgkRingName	= "FMD%c";
const Char_t* AliFMDSimulator::fgkTopHCName	= "F%d%cI";
const Char_t* AliFMDSimulator::fgkBotHCName	= "F%d%cJ";
const Char_t* AliFMDSimulator::fgkTopIHCName	= "F%d%cK";
const Char_t* AliFMDSimulator::fgkBotIHCName	= "F%d%cL";
const Char_t* AliFMDSimulator::fgkNoseName      = "F3SN";
const Char_t* AliFMDSimulator::fgkBackName      = "F3SB";
const Char_t* AliFMDSimulator::fgkBeamName      = "F3SL";
const Char_t* AliFMDSimulator::fgkFlangeName    = "F3SF";

//____________________________________________________________________
AliFMDSimulator::AliFMDSimulator() 
  : fFMD(0), 
    fDetailed(kFALSE),
    fInnerId(-1),
    fOuterId(-1)
{
  // Default constructor
}

//____________________________________________________________________
AliFMDSimulator::AliFMDSimulator(AliFMD* fmd, Bool_t detailed) 
  : TTask("FMDsimulator", "Forward Multiplicity Detector Simulator"), 
    fFMD(fmd), 
    fDetailed(detailed),
    fInnerId(-1),
    fOuterId(-1)
{
  // Normal constructor
  // 
  // Parameters: 
  // 
  //      fmd		Pointer to AliFMD object 
  //      detailed      Whether to make a detailed simulation or not 
  // 
}


//____________________________________________________________________
void
AliFMDSimulator::DefineMaterials() 
{
  // Define the materials and tracking mediums needed by the FMD
  // simulation.   These mediums are made by sending the messages
  // AliMaterial, AliMixture, and AliMedium to the passed AliModule
  // object module.   The defined mediums are 
  // 
  //	FMD Si$		Silicon (active medium in sensors)
  //	FMD C$		Carbon fibre (support cone for FMD3 and vacuum pipe)
  //	FMD Al$		Aluminium (honeycomb support plates)
  //	FMD PCB$	Printed Circuit Board (FEE board with VA1_ALICE)
  //	FMD Chip$	Electronics chips (currently not used)
  //	FMD Air$	Air (Air in the FMD)
  //	FMD Plastic$	Plastic (Support legs for the hybrid cards)
  //
  // Pointers to TGeoMedium objects are retrived from the TGeoManager
  // singleton.  These pointers are later used when setting up the
  // geometry 
  AliDebug(10, "\tCreating materials");
  // Get pointer to geometry singleton object. 
  AliFMDGeometry* geometry = AliFMDGeometry::Instance();
  geometry->Init();
  
  Int_t    id;
  Double_t a                = 0;
  Double_t z                = 0;
  Double_t density          = 0;
  Double_t radiationLength  = 0;
  Double_t absorbtionLength = 999;
  Int_t    fieldType        = gAlice->Field()->Integ();     // Field type 
  Double_t maxField         = gAlice->Field()->Max();     // Field max.
  Double_t maxBending       = 0;     // Max Angle
  Double_t maxStepSize      = 0.001; // Max step size 
  Double_t maxEnergyLoss    = 1;     // Max Delta E
  Double_t precision        = 0.001; // Precision
  Double_t minStepSize      = 0.001; // Minimum step size 
 
  // Silicon 
  a                = 28.0855;
  z                = 14.;
  density          = geometry->GetSiDensity();
  radiationLength  = 9.36;
  maxBending       = 1;
  maxStepSize      = .001;
  precision        = .001;
  minStepSize      = .001;
  id               = kSiId;
  fFMD->AliMaterial(id, "Si$", 
		      a, z, density, radiationLength, absorbtionLength);
  fFMD->AliMedium(kSiId, "Si$",
		    id,1,fieldType,maxField,maxBending,
		    maxStepSize,maxEnergyLoss,precision,minStepSize);
  

  // Carbon 
  a                = 12.011;
  z                = 6.;
  density          = 2.265;
  radiationLength  = 18.8;
  maxBending       = 10;
  maxStepSize      = .01;
  precision        = .003;
  minStepSize      = .003;
  id               = kCarbonId;
  fFMD->AliMaterial(id, "Carbon$", 
		      a, z, density, radiationLength, absorbtionLength);
  fFMD->AliMedium(kCarbonId, "Carbon$",
		    id,0,fieldType,maxField,maxBending,
		    maxStepSize,maxEnergyLoss,precision,minStepSize);

  // Aluminum
  a                = 26.981539;
  z                = 13.;
  density          = 2.7;
  radiationLength  = 8.9;
  id               = kAlId;
  fFMD->AliMaterial(id, "Aluminum$", 
		      a, z, density, radiationLength, absorbtionLength);
  fFMD->AliMedium(kAlId, "Aluminum$", 
		    id, 0, fieldType, maxField, maxBending,
		    maxStepSize, maxEnergyLoss, precision, minStepSize);
  
  
  // Silicon chip 
  {
    Float_t as[] = { 12.0107,      14.0067,      15.9994,
		     1.00794,      28.0855,     107.8682 };
    Float_t zs[] = {  6.,           7.,           8.,
		      1.,          14.,          47. };
    Float_t ws[] = {  0.039730642,  0.001396798,  0.01169634,
		      0.004367771,  0.844665,     0.09814344903 };
    density = 2.36436;
    maxBending       = 10;
    maxStepSize      = .01;
    precision        = .003;
    minStepSize      = .003;
    id = kSiChipId;
    fFMD->AliMixture(id, "Si Chip$", as, zs, density, 6, ws);
    fFMD->AliMedium(kSiChipId, "Si Chip$", 
		      id, 0, fieldType, maxField, maxBending, 
		      maxStepSize, maxEnergyLoss, precision, minStepSize);
  }
  
#if 0
  // Kaption
  {
    Float_t as[] = { 1.00794,  12.0107,  14.010,   15.9994};
    Float_t zs[] = { 1.,        6.,       7.,       8.};
    Float_t ws[] = { 0.026362,  0.69113,  0.07327,  0.209235};
    density          = 1.42;
    maxBending       = 1;
    maxStepSize      = .001;
    precision        = .001;
    minStepSize      = .001;
    id               = KaptionId;
    fFMD->AliMixture(id, "Kaption$", as, zs, density, 4, ws);
    fFMD->AliMedium(kAlId, "Kaption$",
		      id,0,fieldType,maxField,maxBending,
		      maxStepSize,maxEnergyLoss,precision,minStepSize);
  }
#endif

  // Air
  {
    Float_t as[] = { 12.0107, 14.0067,   15.9994,  39.948 };
    Float_t zs[] = {  6.,      7.,       8.,       18. };
    Float_t ws[] = { 0.000124, 0.755267, 0.231781, 0.012827 }; 
    density      = .00120479;
    maxBending   = 1;
    maxStepSize  = .001;
    precision    = .001;
    minStepSize  = .001;
    id           = kAirId;
    fFMD->AliMixture(id, "Air$", as, zs, density, 4, ws);
    fFMD->AliMedium(kAirId, "Air$", 
		      id,0,fieldType,maxField,maxBending,
		      maxStepSize,maxEnergyLoss,precision,minStepSize);
  }
  
  // PCB
  {
    Float_t zs[] = { 14.,         20.,         13.,         12.,
		      5.,         22.,         11.,         19.,
		     26.,          9.,          8.,          6.,
		      7.,          1.};
    Float_t as[] = { 28.0855,     40.078,      26.981538,   24.305, 
		     10.811,      47.867,      22.98977,    39.0983,
		     55.845,      18.9984,     15.9994,     12.0107,
		     14.0067,      1.00794};
    Float_t ws[] = {  0.15144894,  0.08147477,  0.04128158,  0.00904554, 
		      0.01397570,  0.00287685,  0.00445114,  0.00498089,
		      0.00209828,  0.00420000,  0.36043788,  0.27529426,
		      0.01415852,  0.03427566};
    density      = 1.8;
    maxBending   = 1;
    maxStepSize  = .001;
    precision    = .001;
    minStepSize  = .001;
    id           = kPcbId;
    fFMD->AliMixture(id, "PCB$", as, zs, density, 14, ws);
    fFMD->AliMedium(kPcbId, "PCB$", 
		      id,0,fieldType,maxField,maxBending,
		      maxStepSize,maxEnergyLoss,precision,minStepSize);
  }
  
  // Plastic 
  {
    Float_t as[] = { 1.01, 12.01 };
    Float_t zs[] = { 1.,   6.    };
    Float_t ws[] = { 1.,   1.    };
    density      = 1.03;
    maxBending   = 10;
    maxStepSize  = .01;
    precision    = .003;
    minStepSize  = .003;
    id           = kPlasticId;
    fFMD->AliMixture(id, "Plastic$", as, zs, density, -2, ws);
    fFMD->AliMedium(kPlasticId, "Plastic$", 
		      id,0,fieldType,maxField,maxBending,
		      maxStepSize,maxEnergyLoss,precision,minStepSize);
  }
}

//____________________________________________________________________
void
AliFMDSimulator::Exec(Option_t* /* option */) 
{
  // Member function that is executed each time a hit is made in the
  // FMD.  None-charged particles are ignored.   Dead tracks  are
  // ignored. 
  //
  // The procedure is as follows: 
  // 
  //   - IF NOT track is alive THEN RETURN ENDIF
  //   - IF NOT particle is charged THEN RETURN ENDIF
  //   - IF NOT volume name is "STRI" or "STRO" THEN RETURN ENDIF 
  //   - Get strip number (volume copy # minus 1)
  //   - Get phi division number (mother volume copy #)
  //   - Get module number (grand-mother volume copy #)
  //   - section # = 2 * module # + phi division # - 1
  //   - Get ring Id from volume name 
  //   - Get detector # from grand-grand-grand-mother volume name 
  //   - Get pointer to sub-detector object. 
  //   - Get track position 
  //   - IF track is entering volume AND track is inside real shape THEN
  //   -   Reset energy deposited 
  //   -   Get track momentum 
  //   -   Get particle ID # 
  ///  - ENDIF
  //   - IF track is inside volume AND inside real shape THEN 
  ///  -   Update energy deposited 
  //   - ENDIF 
  //   - IF track is inside real shape AND (track is leaving volume,
  //         or it died, or it is stopped  THEN
  //   -   Create a hit 
  //   - ENDIF
  //     
  TVirtualMC* mc = TVirtualMC::GetMC();
  
  if (!mc->IsTrackAlive()) return;
  if (TMath::Abs(mc->TrackCharge()) <= 0) return;
  
  Int_t copy;
  Int_t vol = mc->CurrentVolID(copy);
  if (vol != fInnerId && vol != fOuterId) {
    AliDebug(15, Form("Not an FMD volume %d '%s' (%d or %d)", 
		      vol, mc->CurrentVolName(), fInnerId, fOuterId));
    return;
  }

  // Check that the track is actually within the active area 
  Bool_t entering = mc->IsTrackEntering();
  Bool_t inside   = mc->IsTrackInside();
  Bool_t out      = (mc->IsTrackExiting()|| mc->IsTrackDisappeared()||
		     mc->IsTrackStop());

  // Reset the energy deposition for this track, and update some of
  // our parameters.
  if (entering) {
    AliDebug(15, "Entering active FMD volume");
    fCurrentDeltaE = 0;

    // Get production vertex and momentum of the track 
    mc->TrackMomentum(fCurrentP);
    mc->TrackPosition(fCurrentV);
    fCurrentPdg = mc->IdFromPDG(mc->TrackPid());
  }
  
  // If the track is inside, then update the energy deposition
  if (inside && fCurrentDeltaE >= 0) 
    AliDebug(15, "Inside active FMD volume");
    fCurrentDeltaE += 1000 * mc->Edep();

  // The track exits the volume, or it disappeared in the volume, or
  // the track is stopped because it no longer fulfills the cuts
  // defined, then we create a hit. 
  if (out && fCurrentDeltaE >= 0) {
    AliDebug(15, Form("Leaving active FMD volume %s", mc->CurrentVolPath()));

    Int_t strip = copy - 1;
    Int_t sectordiv;
    mc->CurrentVolOffID(fSectorOff, sectordiv);
    Int_t module;
    mc->CurrentVolOffID(fModuleOff, module);
    Int_t sector = 2 * module + sectordiv;
    Int_t iring;
    mc->CurrentVolOffID(fRingOff, iring);
    Char_t ring = Char_t(iring);
    Int_t detector;
    mc->CurrentVolOffID(fDetectorOff, detector);


    AliFMDGeometry* fmd = AliFMDGeometry::Instance();
    Double_t rz = fmd->GetDetector(detector)->GetRingZ(ring);
    Int_t    n  = fmd->GetDetector(detector)->GetRing(ring)->GetNSectors();
    if (rz < 0) {
      Int_t s = ((n - sector + n / 2) % n) + 1;
      AliDebug(40, Form("Recalculating sector to %d (=%d-%d+%d/2%%%d+1 z=%f)", 
			s, n, sector, n, n, rz));
      sector = s;
    }
    if (sector < 1 || sector > n) {
      Warning("Step", "sector # %d out of range (0-%d)", sector-1, n-1);
      return;
    }
    sector--;
    fCurrentDeltaE += 1000 * mc->Edep();

    AliDebug(20, Form("Processing hit in FMD%d%c[%2d,%3d]: %f", 
		      detector, ring, sector, strip, fCurrentDeltaE));
    
    fFMD->AddHitByFields(gAlice->GetMCApp()->GetCurrentTrackNumber(),
			 UShort_t(detector), ring, UShort_t(sector), 
			 UShort_t(strip),
			 fCurrentV.X(), fCurrentV.Y(), fCurrentV.Z(),
			 fCurrentP.X(), fCurrentP.Y(), fCurrentP.Z(), 
			 fCurrentDeltaE, fCurrentPdg, fCurrentV.T());
    fCurrentDeltaE = -1;
  }
}


//____________________________________________________________________
//
// EOF
//
Commit	Line	Data
1a1fdef7	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	/* $Id$ */
	17
	18	//____________________________________________________________________
	19	//
	20	// Forward Multiplicity Detector based on Silicon wafers. This class
	21	// contains the base procedures for the Forward Multiplicity detector
	22	// Detector consists of 3 sub-detectors FMD1, FMD2, and FMD3, each of
	23	// which has 1 or 2 rings of silicon sensors.
	24	//
	25	// This is the base class for all FMD manager classes.
	26	//
	27	// The actual code is done by various separate classes. Below is
	28	// diagram showing the relationship between the various FMD classes
	29	// that handles the simulation
	30	//
	31	// +--------+ 1 +-----------------+
	32	// \| AliFMD \|<>-----\| AliFMDSimulator \|
	33	// +--------+ +-----------------+
	34	// ^
	35	// \|
	36	// +-------------+-------------+
	37	// \| \|
	38	// +--------------------+ +-------------------+
	39	// \| AliFMDGeoSimulator \| \| AliFMDG3Simulator \|
	40	// +--------------------+ +---------+---------+
	41	//
	42	//
	43	// * AliFMD
	44	// This defines the interface for the various parts of AliROOT that
	45	// uses the FMD, like AliFMDSimulator, AliFMDDigitizer,
	46	// AliFMDReconstructor, and so on.
	47	//
	48	// * AliFMDSimulator
	49	// This is the base class for the FMD simulation tasks. The
	50	// simulator tasks are responsible to implment the geoemtry, and
	51	// process hits.
	52	//
	53	// * AliFMDGeoSimulator
	54	// This is a concrete implementation of the AliFMDSimulator that
	55	// uses the TGeo classes directly only. This defines the active
	56	// volume as an ONLY XTRU shape with a divided MANY TUBS shape
	57	// inside to implement the particular shape of the silicon
	58	// sensors.
	59	//
	60	// * AliFMDG3Simulator
	61	// This is a concrete implementation of the AliFMDSimulator that
	62	// uses the TVirtualMC interface with GEANT 3.21-like messages.
	63	// This implements the active volume as a divided TUBS shape. Hits
	64	// in the corners should be cut away at run time (but currently
65	// isn't).
66	//
67	#include "AliFMDSimulator.h" // ALIFMDSIMULATOR_H
68	#include "AliFMDGeometry.h" // ALIFMDGEOMETRY_H
69	#include "AliFMDDetector.h" // ALIFMDDETECTOR_H
70	#include "AliFMDRing.h" // ALIFMDRING_H
71	#include "AliFMD1.h" // ALIFMD1_H
72	#include "AliFMD2.h" // ALIFMD2_H
73	#include "AliFMD3.h" // ALIFMD3_H
74	#include "AliFMD.h" // ALIFMD_H
75	#include <AliRun.h> // ALIRUN_H
76	#include <AliMC.h> // ALIMC_H
77	#include <AliMagF.h> // ALIMAGF_H
78	#include <AliLog.h> // ALILOG_H
79	#include <TGeoVolume.h> // ROOT_TGeoVolume
80	#include <TGeoTube.h> // ROOT_TGeoTube
81	#include <TGeoPcon.h> // ROOT_TGeoPcon
82	#include <TGeoMaterial.h> // ROOT_TGeoMaterial
83	#include <TGeoMedium.h> // ROOT_TGeoMedium
84	#include <TGeoXtru.h> // ROOT_TGeoXtru
85	#include <TGeoPolygon.h> // ROOT_TGeoPolygon
86	#include <TGeoTube.h> // ROOT_TGeoTube
87	#include <TGeoManager.h> // ROOT_TGeoManager
88	#include <TTree.h> // ROOT_TTree
89	#include <TParticle.h> // ROOT_TParticle
90	#include <TLorentzVector.h> // ROOT_TLorentzVector
91	#include <TVector2.h> // ROOT_TVector2
92	#include <TVector3.h> // ROOT_TVector3
93	#include <TVirtualMC.h> // ROOT_TVirtualMC
94	#include <TArrayD.h> // ROOT_TArrayD
95
96	//====================================================================
97	ClassImp(AliFMDSimulator)
98	#if 0
99	; // This is here to keep Emacs for indenting the next line
100	#endif
101
102	//____________________________________________________________________
103	const Char_t* AliFMDSimulator::fgkActiveName = "F%cAC";
104	const Char_t* AliFMDSimulator::fgkSectorName = "F%cSE";
105	const Char_t* AliFMDSimulator::fgkStripName = "F%cST";
106	const Char_t* AliFMDSimulator::fgkModuleName = "F%cMO";
107	const Char_t* AliFMDSimulator::fgkPCBName = "F%cP%c";
108	const Char_t* AliFMDSimulator::fgkLongLegName = "F%cLL";
109	const Char_t* AliFMDSimulator::fgkShortLegName = "F%cSL";
110	const Char_t* AliFMDSimulator::fgkFrontVName = "F%cFV";
111	const Char_t* AliFMDSimulator::fgkBackVName = "F%cBV";
112	const Char_t* AliFMDSimulator::fgkRingName = "FMD%c";
113	const Char_t* AliFMDSimulator::fgkTopHCName = "F%d%cI";
114	const Char_t* AliFMDSimulator::fgkBotHCName = "F%d%cJ";
115	const Char_t* AliFMDSimulator::fgkTopIHCName = "F%d%cK";
116	const Char_t* AliFMDSimulator::fgkBotIHCName = "F%d%cL";
117	const Char_t* AliFMDSimulator::fgkNoseName = "F3SN";
118	const Char_t* AliFMDSimulator::fgkBackName = "F3SB";
119	const Char_t* AliFMDSimulator::fgkBeamName = "F3SL";
120	const Char_t* AliFMDSimulator::fgkFlangeName = "F3SF";
121
122	//____________________________________________________________________
123	AliFMDSimulator::AliFMDSimulator()
124	: fFMD(0),
125	fDetailed(kFALSE),
126	fInnerId(-1),
127	fOuterId(-1)
128	{
129	// Default constructor
130	}
131
132	//____________________________________________________________________
133	AliFMDSimulator::AliFMDSimulator(AliFMD* fmd, Bool_t detailed)
134	: TTask("FMDsimulator", "Forward Multiplicity Detector Simulator"),
135	fFMD(fmd),
136	fDetailed(detailed),
137	fInnerId(-1),
138	fOuterId(-1)
139	{
140	// Normal constructor
141	//
142	// Parameters:
143	//
144	// fmd Pointer to AliFMD object
145	// detailed Whether to make a detailed simulation or not
146	//
147	}
148
149
150	//____________________________________________________________________
151	void
152	AliFMDSimulator::DefineMaterials()
153	{
154	// Define the materials and tracking mediums needed by the FMD
155	// simulation. These mediums are made by sending the messages
156	// AliMaterial, AliMixture, and AliMedium to the passed AliModule
157	// object module. The defined mediums are
158	//
159	// FMD Si$ Silicon (active medium in sensors)
160	// FMD C$ Carbon fibre (support cone for FMD3 and vacuum pipe)
161	// FMD Al$ Aluminium (honeycomb support plates)
162	// FMD PCB$ Printed Circuit Board (FEE board with VA1_ALICE)
163	// FMD Chip$ Electronics chips (currently not used)
164	// FMD Air$ Air (Air in the FMD)
165	// FMD Plastic$ Plastic (Support legs for the hybrid cards)
166	//
167	// Pointers to TGeoMedium objects are retrived from the TGeoManager
168	// singleton. These pointers are later used when setting up the
169	// geometry
170	AliDebug(10, "\tCreating materials");
171	// Get pointer to geometry singleton object.
172	AliFMDGeometry* geometry = AliFMDGeometry::Instance();
173	geometry->Init();
174
175	Int_t id;
176	Double_t a = 0;
177	Double_t z = 0;
178	Double_t density = 0;
179	Double_t radiationLength = 0;
180	Double_t absorbtionLength = 999;
181	Int_t fieldType = gAlice->Field()->Integ(); // Field type
182	Double_t maxField = gAlice->Field()->Max(); // Field max.
183	Double_t maxBending = 0; // Max Angle
184	Double_t maxStepSize = 0.001; // Max step size
185	Double_t maxEnergyLoss = 1; // Max Delta E
186	Double_t precision = 0.001; // Precision
187	Double_t minStepSize = 0.001; // Minimum step size
188
189	// Silicon
190	a = 28.0855;
191	z = 14.;
192	density = geometry->GetSiDensity();
193	radiationLength = 9.36;
194	maxBending = 1;
195	maxStepSize = .001;
196	precision = .001;
197	minStepSize = .001;
198	id = kSiId;
4a9de4af	199	fFMD->AliMaterial(id, "Si$",
1a1fdef7	200	a, z, density, radiationLength, absorbtionLength);
4a9de4af	201	fFMD->AliMedium(kSiId, "Si$",
1a1fdef7	202	id,1,fieldType,maxField,maxBending,
	203	maxStepSize,maxEnergyLoss,precision,minStepSize);
	204
	205
	206	// Carbon
	207	a = 12.011;
	208	z = 6.;
	209	density = 2.265;
	210	radiationLength = 18.8;
	211	maxBending = 10;
	212	maxStepSize = .01;
	213	precision = .003;
	214	minStepSize = .003;
	215	id = kCarbonId;
4a9de4af	216	fFMD->AliMaterial(id, "Carbon$",
1a1fdef7	217	a, z, density, radiationLength, absorbtionLength);
4a9de4af	218	fFMD->AliMedium(kCarbonId, "Carbon$",
1a1fdef7	219	id,0,fieldType,maxField,maxBending,
	220	maxStepSize,maxEnergyLoss,precision,minStepSize);
	221
	222	// Aluminum
	223	a = 26.981539;
	224	z = 13.;
	225	density = 2.7;
	226	radiationLength = 8.9;
	227	id = kAlId;
4a9de4af	228	fFMD->AliMaterial(id, "Aluminum$",
1a1fdef7	229	a, z, density, radiationLength, absorbtionLength);
4a9de4af	230	fFMD->AliMedium(kAlId, "Aluminum$",
1a1fdef7	231	id, 0, fieldType, maxField, maxBending,
	232	maxStepSize, maxEnergyLoss, precision, minStepSize);
	233
	234
	235	// Silicon chip
	236	{
	237	Float_t as[] = { 12.0107, 14.0067, 15.9994,
	238	1.00794, 28.0855, 107.8682 };
	239	Float_t zs[] = { 6., 7., 8.,
	240	1., 14., 47. };
	241	Float_t ws[] = { 0.039730642, 0.001396798, 0.01169634,
	242	0.004367771, 0.844665, 0.09814344903 };
	243	density = 2.36436;
	244	maxBending = 10;
	245	maxStepSize = .01;
	246	precision = .003;
	247	minStepSize = .003;
	248	id = kSiChipId;
4a9de4af	249	fFMD->AliMixture(id, "Si Chip$", as, zs, density, 6, ws);
4a9de4af	250	fFMD->AliMedium(kSiChipId, "Si Chip$",
1a1fdef7	251	id, 0, fieldType, maxField, maxBending,
	252	maxStepSize, maxEnergyLoss, precision, minStepSize);
	253	}
	254
	255	#if 0
	256	// Kaption
	257	{
	258	Float_t as[] = { 1.00794, 12.0107, 14.010, 15.9994};
	259	Float_t zs[] = { 1., 6., 7., 8.};
	260	Float_t ws[] = { 0.026362, 0.69113, 0.07327, 0.209235};
	261	density = 1.42;
	262	maxBending = 1;
	263	maxStepSize = .001;
	264	precision = .001;
	265	minStepSize = .001;
	266	id = KaptionId;
4a9de4af	267	fFMD->AliMixture(id, "Kaption$", as, zs, density, 4, ws);
4a9de4af	268	fFMD->AliMedium(kAlId, "Kaption$",
1a1fdef7	269	id,0,fieldType,maxField,maxBending,
	270	maxStepSize,maxEnergyLoss,precision,minStepSize);
	271	}
	272	#endif
	273
	274	// Air
	275	{
	276	Float_t as[] = { 12.0107, 14.0067, 15.9994, 39.948 };
	277	Float_t zs[] = { 6., 7., 8., 18. };
	278	Float_t ws[] = { 0.000124, 0.755267, 0.231781, 0.012827 };
	279	density = .00120479;
	280	maxBending = 1;
	281	maxStepSize = .001;
	282	precision = .001;
	283	minStepSize = .001;
	284	id = kAirId;
4a9de4af	285	fFMD->AliMixture(id, "Air$", as, zs, density, 4, ws);
4a9de4af	286	fFMD->AliMedium(kAirId, "Air$",
1a1fdef7	287	id,0,fieldType,maxField,maxBending,
	288	maxStepSize,maxEnergyLoss,precision,minStepSize);
	289	}
	290
	291	// PCB
	292	{
	293	Float_t zs[] = { 14., 20., 13., 12.,
	294	5., 22., 11., 19.,
	295	26., 9., 8., 6.,
	296	7., 1.};
	297	Float_t as[] = { 28.0855, 40.078, 26.981538, 24.305,
	298	10.811, 47.867, 22.98977, 39.0983,
	299	55.845, 18.9984, 15.9994, 12.0107,
	300	14.0067, 1.00794};
	301	Float_t ws[] = { 0.15144894, 0.08147477, 0.04128158, 0.00904554,
	302	0.01397570, 0.00287685, 0.00445114, 0.00498089,
	303	0.00209828, 0.00420000, 0.36043788, 0.27529426,
	304	0.01415852, 0.03427566};
	305	density = 1.8;
	306	maxBending = 1;
	307	maxStepSize = .001;
	308	precision = .001;
	309	minStepSize = .001;
	310	id = kPcbId;
4a9de4af	311	fFMD->AliMixture(id, "PCB$", as, zs, density, 14, ws);
4a9de4af	312	fFMD->AliMedium(kPcbId, "PCB$",
1a1fdef7	313	id,0,fieldType,maxField,maxBending,
	314	maxStepSize,maxEnergyLoss,precision,minStepSize);
	315	}
	316
	317	// Plastic
	318	{
	319	Float_t as[] = { 1.01, 12.01 };
	320	Float_t zs[] = { 1., 6. };
	321	Float_t ws[] = { 1., 1. };
	322	density = 1.03;
	323	maxBending = 10;
	324	maxStepSize = .01;
	325	precision = .003;
	326	minStepSize = .003;
	327	id = kPlasticId;
4a9de4af	328	fFMD->AliMixture(id, "Plastic$", as, zs, density, -2, ws);
4a9de4af	329	fFMD->AliMedium(kPlasticId, "Plastic$",
1a1fdef7	330	id,0,fieldType,maxField,maxBending,
	331	maxStepSize,maxEnergyLoss,precision,minStepSize);
	332	}
	333	}
	334
	335	//____________________________________________________________________
	336	void
	337	AliFMDSimulator::Exec(Option_t* /* option */)
	338	{
	339	// Member function that is executed each time a hit is made in the
	340	// FMD. None-charged particles are ignored. Dead tracks are
	341	// ignored.
	342	//
	343	// The procedure is as follows:
	344	//
	345	// - IF NOT track is alive THEN RETURN ENDIF
	346	// - IF NOT particle is charged THEN RETURN ENDIF
	347	// - IF NOT volume name is "STRI" or "STRO" THEN RETURN ENDIF
	348	// - Get strip number (volume copy # minus 1)
	349	// - Get phi division number (mother volume copy #)
	350	// - Get module number (grand-mother volume copy #)
	351	// - section # = 2 * module # + phi division # - 1
	352	// - Get ring Id from volume name
	353	// - Get detector # from grand-grand-grand-mother volume name
	354	// - Get pointer to sub-detector object.
	355	// - Get track position
	356	// - IF track is entering volume AND track is inside real shape THEN
	357	// - Reset energy deposited
	358	// - Get track momentum
	359	// - Get particle ID #
	360	/// - ENDIF
	361	// - IF track is inside volume AND inside real shape THEN
	362	/// - Update energy deposited
	363	// - ENDIF
	364	// - IF track is inside real shape AND (track is leaving volume,
	365	// or it died, or it is stopped THEN
	366	// - Create a hit
	367	// - ENDIF
	368	//
	369	TVirtualMC* mc = TVirtualMC::GetMC();
	370
	371	if (!mc->IsTrackAlive()) return;
	372	if (TMath::Abs(mc->TrackCharge()) <= 0) return;
	373
	374	Int_t copy;
	375	Int_t vol = mc->CurrentVolID(copy);
	376	if (vol != fInnerId && vol != fOuterId) {
	377	AliDebug(15, Form("Not an FMD volume %d '%s' (%d or %d)",
	378	vol, mc->CurrentVolName(), fInnerId, fOuterId));
	379	return;
	380	}
	381
	382	// Check that the track is actually within the active area
	383	Bool_t entering = mc->IsTrackEntering();
	384	Bool_t inside = mc->IsTrackInside();
	385	Bool_t out = (mc->IsTrackExiting()\|\| mc->IsTrackDisappeared()\|\|
	386	mc->IsTrackStop());
	387
	388	// Reset the energy deposition for this track, and update some of
	389	// our parameters.
	390	if (entering) {
	391	AliDebug(15, "Entering active FMD volume");
	392	fCurrentDeltaE = 0;
	393
394	// Get production vertex and momentum of the track
395	mc->TrackMomentum(fCurrentP);
396	mc->TrackPosition(fCurrentV);
397	fCurrentPdg = mc->IdFromPDG(mc->TrackPid());
398	}
399
400	// If the track is inside, then update the energy deposition
401	if (inside && fCurrentDeltaE >= 0)
402	AliDebug(15, "Inside active FMD volume");
403	fCurrentDeltaE += 1000 * mc->Edep();
404
405	// The track exits the volume, or it disappeared in the volume, or
406	// the track is stopped because it no longer fulfills the cuts
407	// defined, then we create a hit.
408	if (out && fCurrentDeltaE >= 0) {
409	AliDebug(15, Form("Leaving active FMD volume %s", mc->CurrentVolPath()));
410
411	Int_t strip = copy - 1;
412	Int_t sectordiv;
413	mc->CurrentVolOffID(fSectorOff, sectordiv);
414	Int_t module;
415	mc->CurrentVolOffID(fModuleOff, module);
416	Int_t sector = 2 * module + sectordiv;
417	Int_t iring;
418	mc->CurrentVolOffID(fRingOff, iring);
419	Char_t ring = Char_t(iring);
420	Int_t detector;
421	mc->CurrentVolOffID(fDetectorOff, detector);
422
423
424	AliFMDGeometry* fmd = AliFMDGeometry::Instance();
425	Double_t rz = fmd->GetDetector(detector)->GetRingZ(ring);
426	Int_t n = fmd->GetDetector(detector)->GetRing(ring)->GetNSectors();
427	if (rz < 0) {
428	Int_t s = ((n - sector + n / 2) % n) + 1;
429	AliDebug(40, Form("Recalculating sector to %d (=%d-%d+%d/2%%%d+1 z=%f)",
430	s, n, sector, n, n, rz));
431	sector = s;
432	}
433	if (sector < 1 \|\| sector > n) {
434	Warning("Step", "sector # %d out of range (0-%d)", sector-1, n-1);
435	return;
436	}
437	sector--;
438	fCurrentDeltaE += 1000 * mc->Edep();
439
440	AliDebug(20, Form("Processing hit in FMD%d%c[%2d,%3d]: %f",
441	detector, ring, sector, strip, fCurrentDeltaE));
442
69b696b9	443	fFMD->AddHitByFields(gAlice->GetMCApp()->GetCurrentTrackNumber(),
	444	UShort_t(detector), ring, UShort_t(sector),
	445	UShort_t(strip),
	446	fCurrentV.X(), fCurrentV.Y(), fCurrentV.Z(),
	447	fCurrentP.X(), fCurrentP.Y(), fCurrentP.Z(),
	448	fCurrentDeltaE, fCurrentPdg, fCurrentV.T());
1a1fdef7	449	fCurrentDeltaE = -1;
	450	}
	451	}
	452
	453
	454
	455	//____________________________________________________________________
	456	//
	457	// EOF
	458	//