Bringing the OB geometry up to date. Macros to produce material budget plots. General...
authorbelikov <Iouri.Belikov@cern.ch>
Wed, 26 Nov 2014 14:52:57 +0000 (15:52 +0100)
committerbelikov <Iouri.Belikov@cern.ch>
Wed, 26 Nov 2014 14:52:57 +0000 (15:52 +0100)
ITS/UPGRADE/AliITSUv1.cxx
ITS/UPGRADE/AliITSUv1.h
ITS/UPGRADE/AliITSUv1Layer.cxx
ITS/UPGRADE/AliITSUv1Layer.h
ITS/UPGRADE/testITSUv1/CreateITSUv1.C
ITS/UPGRADE/testITSUv1/CreateITSUv1.C_template [new file with mode: 0644]
ITS/UPGRADE/testITSUv1/GetMaterialBudget.C [new file with mode: 0644]
ITS/UPGRADE/testITSUv1/MakeMatBudPlots.C [new file with mode: 0644]
ITS/UPGRADE/testITSUv1/MaterialBudget_README [new file with mode: 0644]
ITS/UPGRADE/testITSUv1/runMatBud.sh [new file with mode: 0755]

index a5a715b..d1f4332 100644 (file)
@@ -72,10 +72,10 @@ AliITSUv1::AliITSUv1()
   ,fLayZLength(0)
   ,fStavPerLay(0)
   ,fUnitPerStave(0)
-  ,fStaveThick(0)
+  ,fChipThick(0)
   ,fStaveWidth(0)
   ,fStaveTilt(0)
-  ,fDetThick(0)
+  ,fSensThick(0)
   ,fChipTypeID(0)
   ,fBuildLevel(0)
   ,fUpGeom(0)
@@ -105,10 +105,10 @@ AliITSUv1::AliITSUv1(const char *title, Int_t nlay)
   ,fLayZLength(0)
   ,fStavPerLay(0)
   ,fUnitPerStave(0)
-  ,fStaveThick(0)
+  ,fChipThick(0)
   ,fStaveWidth(0)
   ,fStaveTilt(0)
-  ,fDetThick(0)
+  ,fSensThick(0)
   ,fChipTypeID(0)
   ,fBuildLevel(0)
   ,fUpGeom(0)
@@ -132,10 +132,10 @@ AliITSUv1::AliITSUv1(const char *title, Int_t nlay)
   fLayZLength   = new Double_t[fNLayers];
   fStavPerLay   = new Int_t[fNLayers];
   fUnitPerStave = new Int_t[fNLayers];
-  fStaveThick   = new Double_t[fNLayers];
+  fChipThick    = new Double_t[fNLayers];
   fStaveWidth   = new Double_t[fNLayers];
   fStaveTilt    = new Double_t[fNLayers];
-  fDetThick     = new Double_t[fNLayers];
+  fSensThick    = new Double_t[fNLayers];
   fChipTypeID   = new UInt_t[fNLayers];
   fBuildLevel   = new Int_t[fNLayers];
 
@@ -150,7 +150,7 @@ AliITSUv1::AliITSUv1(const char *title, Int_t nlay)
       fStavPerLay[j]   = 0;
       fUnitPerStave[j] = 0;
       fStaveWidth[j]   = 0.;
-      fDetThick[j]     = 0.;
+      fSensThick[j]    = 0.;
       fChipTypeID[j]   = 0;
       fBuildLevel[j]   = 0;
       fUpGeom[j]       = 0;
@@ -173,10 +173,10 @@ AliITSUv1::~AliITSUv1() {
   delete [] fLayZLength;
   delete [] fStavPerLay;
   delete [] fUnitPerStave;
-  delete [] fStaveThick;
+  delete [] fChipThick;
   delete [] fStaveWidth;
   delete [] fStaveTilt;
-  delete [] fDetThick;
+  delete [] fSensThick;
   delete [] fChipTypeID;
   delete [] fBuildLevel;
   delete [] fUpGeom;
@@ -355,17 +355,17 @@ void AliITSUv1::CreateGeometry() {
     if (fLayZLength[j] <= 0)                 AliFatal(Form("Wrong layer length for layer %d (%f)",j,fLayZLength[j]));
     if (fStavPerLay[j] <= 0)                 AliFatal(Form("Wrong number of staves for layer %d (%d)",j,fStavPerLay[j]));
     if (fUnitPerStave[j] <= 0)               AliFatal(Form("Wrong number of chips for layer %d (%d)",j,fUnitPerStave[j]));
-    if (fStaveThick[j] < 0)                  AliFatal(Form("Wrong stave thickness for layer %d (%f)",j,fStaveThick[j]));
+    if (fChipThick[j] < 0)                   AliFatal(Form("Wrong chip thickness for layer %d (%f)",j,fChipThick[j]));
     if (fLayTurbo[j] && fStaveWidth[j] <= 0) AliFatal(Form("Wrong stave width for layer %d (%f)",j,fStaveWidth[j]));
-    if (fDetThick[j] < 0)                    AliFatal(Form("Wrong chip thickness for layer %d (%f)",j,fDetThick[j]));
+    if (fSensThick[j] < 0)                   AliFatal(Form("Wrong sensor thickness for layer %d (%f)",j,fSensThick[j]));
     //
     if (j > 0) {
       if (fLayRadii[j]<=fLayRadii[j-1])      AliFatal(Form("Layer %d radius (%f) is smaller than layer %d radius (%f)",
                           j,fLayRadii[j],j-1,fLayRadii[j-1]));
     } // if (j > 0)
 
-    if (fStaveThick[j] == 0) AliInfo(Form("Stave thickness for layer %d not set, using default",j));
-    if (fDetThick[j]   == 0) AliInfo(Form("Chip thickness for layer %d not set, using default",j));
+    if (fChipThick[j] == 0) AliInfo(Form("Chip thickness for layer %d not set, using default",j));
+    if (fSensThick[j] == 0) AliInfo(Form("Sensor thickness for layer %d not set, using default",j));
 
   } // for (Int_t j=0; j<fNLayers; j++)
 
@@ -406,8 +406,8 @@ void AliITSUv1::CreateGeometry() {
       fUpGeom[j]->SetStaveModel(fStaveModelOB);
     AliDebug(1,Form("fBuildLevel: %d\n",fBuildLevel[j]));
     //
-    if (fStaveThick[j] != 0) fUpGeom[j]->SetStaveThick(fStaveThick[j]);
-    if (fDetThick[j]   != 0) fUpGeom[j]->SetSensorThick(fDetThick[j]);
+    if (fChipThick[j] != 0) fUpGeom[j]->SetChipThick(fChipThick[j]);
+    if (fSensThick[j] != 0) fUpGeom[j]->SetSensorThick(fSensThick[j]);
     //
     for (int iw=0;iw<fNWrapVol;iw++) {
       if (fLayRadii[j]>fWrapRMin[iw] && fLayRadii[j]<fWrapRMax[iw]) {
@@ -609,8 +609,8 @@ void AliITSUv1::DefineLayer(Int_t nlay, double phi0, Double_t r,
   //          nstav   number of staves
   //          nunit   IB: number of chips per stave
   //                  OB: number of modules per half stave
-  //          lthick  stave thickness (if omitted, defaults to 0)
-  //          dthick  detector thickness (if omitted, defaults to 0)
+  //          lthick  chip thickness (if omitted, defaults to 0)
+  //          dthick  sensor thickness (if omitted, defaults to 0)
   //          dettypeID  ??
   //          buildLevel (if 0, all geometry is build, used for material budget studies)
   // Outputs:
@@ -632,8 +632,8 @@ void AliITSUv1::DefineLayer(Int_t nlay, double phi0, Double_t r,
   fLayZLength[nlay] = zlen;
   fStavPerLay[nlay] = nstav;
   fUnitPerStave[nlay] = nunit;
-  fStaveThick[nlay] = lthick;
-  fDetThick[nlay] = dthick;
+  fChipThick[nlay] = lthick;
+  fSensThick[nlay] = dthick;
   fChipTypeID[nlay] = dettypeID;
   fBuildLevel[nlay] = buildLevel;
     
@@ -657,8 +657,8 @@ void AliITSUv1::DefineLayerTurbo(Int_t nlay, Double_t phi0, Double_t r, Double_t
   //                  OB: number of modules per half stave
   //          width   stave width
   //          tilt    layer tilt angle (degrees)
-  //          lthick  stave thickness (if omitted, defaults to 0)
-  //          dthick  detector thickness (if omitted, defaults to 0)
+  //          lthick  chip thickness (if omitted, defaults to 0)
+  //          dthick  sensor thickness (if omitted, defaults to 0)
   //          dettypeID  ??
   //          buildLevel (if 0, all geometry is build, used for material budget studies)
   // Outputs:
@@ -680,10 +680,10 @@ void AliITSUv1::DefineLayerTurbo(Int_t nlay, Double_t phi0, Double_t r, Double_t
   fLayZLength[nlay] = zlen;
   fStavPerLay[nlay] = nstav;
   fUnitPerStave[nlay] = nunit;
-  fStaveThick[nlay] = lthick;
+  fChipThick[nlay] = lthick;
   fStaveWidth[nlay] = width;
   fStaveTilt[nlay] = tilt;
-  fDetThick[nlay] = dthick;
+  fSensThick[nlay] = dthick;
   fChipTypeID[nlay] = dettypeID;
   fBuildLevel[nlay] = buildLevel;
 
@@ -709,8 +709,8 @@ void AliITSUv1::GetLayerParameters(Int_t nlay, Double_t &phi0,
   //                  OB: number of modules per half stave
   //          width   stave width
   //          tilt    stave tilt angle
-  //          lthick  stave thickness
-  //          dthick  detector thickness
+  //          lthick  chip thickness
+  //          dthick  sensor thickness
   //          dettype detector type
   // Return:
   //   none.
@@ -727,8 +727,8 @@ void AliITSUv1::GetLayerParameters(Int_t nlay, Double_t &phi0,
   nmod   = fUnitPerStave[nlay];
   width  = fStaveWidth[nlay];
   tilt   = fStaveTilt[nlay];
-  lthick = fStaveThick[nlay];
-  dthick = fDetThick[nlay];
+  lthick = fChipThick[nlay];
+  dthick = fSensThick[nlay];
   dettype= fChipTypeID[nlay];
 }
 
index 9cc8fb7..024a827 100644 (file)
@@ -35,7 +35,8 @@ class AliITSUv1 : public AliITSU {
     kIBModel3=5,
     kOBModelDummy=6,
     kOBModel0=7,
-    kOBModel1=8 
+    kOBModel1=8, 
+    kOBModel2=9 
   } AliITSUModel_t;
   
 
@@ -94,10 +95,10 @@ class AliITSUv1 : public AliITSU {
   Double_t *fLayZLength;     // Vector of layer length along Z
   Int_t    *fStavPerLay;     // Vector of number of staves per layer
   Int_t    *fUnitPerStave;   // Vector of number of "units" per stave
-  Double_t *fStaveThick;     // Vector of stave thicknesses
+  Double_t *fChipThick;      // Vector of chip thicknesses
   Double_t *fStaveWidth;     // Vector of stave width (only used for turbo)
   Double_t *fStaveTilt;      // Vector of stave tilt (only used for turbo)
-  Double_t *fDetThick;       // Vector of detector thicknesses
+  Double_t *fSensThick;      // Vector of sensor thicknesses
   UInt_t   *fChipTypeID;     // Vector of detector type id
   Int_t    *fBuildLevel;     // Vector of Material Budget Studies
   //  
@@ -107,7 +108,7 @@ class AliITSUv1 : public AliITSU {
   
   // Parameters for the Upgrade geometry
   
-  ClassDef(AliITSUv1,0)                          
+  ClassDef(AliITSUv1,0)
 };
  
 #endif
index 2ee144e..2446174 100644 (file)
@@ -47,8 +47,8 @@ using namespace TMath;
 // General Parameters
 const Int_t    AliITSUv1Layer::fgkNumberOfInnerLayers =   3;
 
-const Double_t AliITSUv1Layer::fgkDefaultSensorThick  = 300*fgkmicron;
-const Double_t AliITSUv1Layer::fgkDefaultStaveThick   =   1*fgkcm;
+const Double_t AliITSUv1Layer::fgkDefaultSensorThick  =  18*fgkmicron;
+const Double_t AliITSUv1Layer::fgkDefaultChipThick    =  50*fgkmicron;
 
 // Inner Barrel Parameters
 const Int_t    AliITSUv1Layer::fgkIBChipsPerRow       =   9;
@@ -64,17 +64,21 @@ const Double_t AliITSUv1Layer::fgkOBModuleGap         =   0.01 *fgkcm;
 const Double_t AliITSUv1Layer::fgkOBChipXGap          =   0.01 *fgkcm;
 const Double_t AliITSUv1Layer::fgkOBChipZGap          =   0.01 *fgkcm;
 const Double_t AliITSUv1Layer::fgkOBFlexCableAlThick  =   0.005*fgkcm;
+const Double_t AliITSUv1Layer::fgkOBFlexCableCuThick  =   0.004*fgkcm;
 const Double_t AliITSUv1Layer::fgkOBFlexCableKapThick =   0.01 *fgkcm;
 const Double_t AliITSUv1Layer::fgkOBBusCableAlThick   =   0.02 *fgkcm;
 const Double_t AliITSUv1Layer::fgkOBBusCableKapThick  =   0.02 *fgkcm;
 const Double_t AliITSUv1Layer::fgkOBColdPlateThick    =   0.012*fgkcm;
 const Double_t AliITSUv1Layer::fgkOBCarbonPlateThick  =   0.012*fgkcm;
-const Double_t AliITSUv1Layer::fgkOBGlueThick         =   0.03 *fgkcm;
+const Double_t AliITSUv1Layer::fgkOBGlueThickM1       =   0.03 *fgkcm;
+const Double_t AliITSUv1Layer::fgkOBGlueThick         =   0.01 *fgkcm;
 const Double_t AliITSUv1Layer::fgkOBModuleZLength     =  21.06 *fgkcm;
 const Double_t AliITSUv1Layer::fgkOBHalfStaveYTrans   =   1.76 *fgkmm;
 const Double_t AliITSUv1Layer::fgkOBHalfStaveXOverlap =   4.3  *fgkmm;
 const Double_t AliITSUv1Layer::fgkOBGraphiteFoilThick =  30.0  *fgkmicron;
-const Double_t AliITSUv1Layer::fgkOBCoolTubeInnerD    =   2.052*fgkmm;
+const Double_t AliITSUv1Layer::fgkOBCarbonFleeceThick =  20.0  *fgkmicron;
+const Double_t AliITSUv1Layer::fgkOBCoolTubeInnerDM1  =   2.052*fgkmm;
+const Double_t AliITSUv1Layer::fgkOBCoolTubeInnerD    =   2.05 *fgkmm;
 const Double_t AliITSUv1Layer::fgkOBCoolTubeThick     =  32.0  *fgkmicron;
 const Double_t AliITSUv1Layer::fgkOBCoolTubeXDist     =  11.1  *fgkmm;
 
@@ -102,7 +106,7 @@ AliITSUv1Layer::AliITSUv1Layer():
   fLayRadius(0),
   fZLength(0),
   fSensorThick(0),
-  fStaveThick(0),
+  fChipThick(0),
   fStaveWidth(0),
   fStaveTilt(0),
   fNStaves(0),
@@ -127,7 +131,7 @@ AliITSUv1Layer::AliITSUv1Layer(Int_t debug):
   fLayRadius(0),
   fZLength(0),
   fSensorThick(0),
-  fStaveThick(0),
+  fChipThick(0),
   fStaveWidth(0),
   fStaveTilt(0),
   fNStaves(0),
@@ -152,7 +156,7 @@ AliITSUv1Layer::AliITSUv1Layer(Int_t lay, Int_t debug):
   fLayRadius(0),
   fZLength(0),
   fSensorThick(0),
-  fStaveThick(0),
+  fChipThick(0),
   fStaveWidth(0),
   fStaveTilt(0),
   fNStaves(0),
@@ -177,7 +181,7 @@ AliITSUv1Layer::AliITSUv1Layer(Int_t lay, Bool_t turbo, Int_t debug):
   fLayRadius(0),
   fZLength(0),
   fSensorThick(0),
-  fStaveThick(0),
+  fChipThick(0),
   fStaveWidth(0),
   fStaveTilt(0),
   fNStaves(0),
@@ -203,7 +207,7 @@ AliITSUv1Layer::AliITSUv1Layer(const AliITSUv1Layer &s):
   fLayRadius(s.fLayRadius),
   fZLength(s.fZLength),
   fSensorThick(s.fSensorThick),
-  fStaveThick(s.fStaveThick),
+  fChipThick(s.fChipThick),
   fStaveWidth(s.fStaveWidth),
   fStaveTilt(s.fStaveTilt),
   fNStaves(s.fNStaves),
@@ -233,7 +237,7 @@ AliITSUv1Layer& AliITSUv1Layer::operator=(const AliITSUv1Layer &s)
   fLayRadius   = s.fLayRadius;
   fZLength     = s.fZLength;
   fSensorThick = s.fSensorThick;
-  fStaveThick  = s.fStaveThick;
+  fChipThick   = s.fChipThick;
   fStaveWidth  = s.fStaveWidth;
   fStaveTilt   = s.fStaveTilt;
   fNStaves     = s.fNStaves;
@@ -286,10 +290,10 @@ void AliITSUv1Layer::CreateLayer(TGeoVolume *moth){
   if (fLayerNumber >= fgkNumberOfInnerLayers && fNModules <= 0)
     AliFatal(Form("Wrong number of modules (%d)",fNModules));
 
-  if (fStaveThick <= 0) {
-    AliInfo(Form("Stave thickness wrong or not set (%f), using default (%f)",
-                fStaveThick,fgkDefaultStaveThick));
-    fStaveThick = fgkDefaultStaveThick;
+  if (fChipThick <= 0) {
+    AliInfo(Form("Chip thickness wrong or not set (%f), using default (%f)",
+                fChipThick,fgkDefaultChipThick));
+    fChipThick = fgkDefaultChipThick;
   }
 
   if (fSensorThick <= 0) {
@@ -298,10 +302,10 @@ void AliITSUv1Layer::CreateLayer(TGeoVolume *moth){
     fSensorThick = fgkDefaultSensorThick;
   }
 
-  if (fSensorThick > fStaveThick) {
-    AliWarning(Form("Sensor thickness (%f) is greater than stave thickness (%f), fixing",
-                fSensorThick,fStaveThick));
-    fSensorThick = fStaveThick;
+  if (fSensorThick > fChipThick) {
+    AliWarning(Form("Sensor thickness (%f) is greater than chip thickness (%f), fixing",
+                fSensorThick,fChipThick));
+    fSensorThick = fChipThick;
   }
 
 
@@ -319,6 +323,8 @@ void AliITSUv1Layer::CreateLayer(TGeoVolume *moth){
 
   snprintf(volname, 30, "%s%d", AliITSUGeomTGeo::GetITSLayerPattern(),fLayerNumber);
   TGeoVolume *layVol = new TGeoVolumeAssembly(volname);
+  layVol->SetTitle(
+          Form("Model %d Build level %d",(Int_t)fStaveModel,fBuildLevel));
   layVol->SetUniqueID(fChipTypeID);
 
 //  layVol->SetVisibility(kFALSE);
@@ -386,6 +392,8 @@ void AliITSUv1Layer::CreateLayerTurbo(TGeoVolume *moth){
 
   snprintf(volname, 30, "%s%d", AliITSUGeomTGeo::GetITSLayerPattern(), fLayerNumber);
   TGeoVolume *layVol = new TGeoVolumeAssembly(volname);
+  layVol->SetTitle(
+          Form("Model %d Build level %d",(Int_t)fStaveModel,fBuildLevel));
   layVol->SetUniqueID(fChipTypeID);
   layVol->SetVisibility(kTRUE);
   layVol->SetLineColor(1);
@@ -441,7 +449,7 @@ TGeoVolume* AliITSUv1Layer::CreateStave(const TGeoManager * /*mgr*/){
   // The stave
   xlen = fLayRadius*Tan(alpha);
   if (fIsTurbo) xlen = 0.5*fStaveWidth;
-  ylen = 0.5*fStaveThick;
+  ylen = 0.5*fChipThick;
   zlen = 0.5*fZLength;
 
   Double_t yplus = 0.46;
@@ -1930,7 +1938,8 @@ TGeoVolume* AliITSUv1Layer::CreateStaveOuterB(const TGeoManager *mgr){
       mechStavVol = CreateStaveModelOuterB0(mgr);
       break;
     case AliITSUv1::kOBModel1:
-      mechStavVol = CreateStaveModelOuterB1(mgr);
+    case AliITSUv1::kOBModel2:
+      mechStavVol = CreateStaveModelOuterB12(mgr);
       break;
     default:
       AliFatal(Form("Unknown stave model %d",fStaveModel));
@@ -1988,7 +1997,7 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB0(const TGeoManager *mgr){
 
   // The chip
   xlen = fgkOBHalfStaveWidth;
-  ylen = 0.5*fStaveThick; // TO BE CHECKED
+  ylen = 0.5*fChipThick;
   zlen = fgkOBModuleZLength/2;
 
   TGeoVolume *chipVol = CreateChipInnerB(xlen, ylen, zlen);
@@ -2032,7 +2041,7 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB0(const TGeoManager *mgr){
 }
 
 //________________________________________________________________________
-TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
+TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB12(const TGeoManager *mgr){
 //
 // Create the mechanical half stave structure
 // for the Outer Barrel as in TDR
@@ -2047,26 +2056,31 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
 // Created:      20 Nov 2013  Anastasia Barbano
 // Updated:      16 Jan 2014  Mario Sitta
 // Updated:      24 Feb 2014  Mario Sitta
+// Updated:      11 Nov 2014  Mario Sitta  Model2
 //
 
 
   // Local parameters
   Double_t yFlex1      = fgkOBFlexCableAlThick;
   Double_t yFlex2      = fgkOBFlexCableKapThick;
-  Double_t flexOverlap = 5; // to be checked
+  Double_t flexOverlap = 5; // to be checked - unused for the time being
   Double_t xHalfSt     = fgkOBHalfStaveWidth/2;
-  Double_t rCoolMin    = fgkOBCoolTubeInnerD/2;
-  Double_t rCoolMax    = rCoolMin + fgkOBCoolTubeThick;
-  Double_t kLay1       = 0.004; // to be checked
+  Double_t kLay1       = fgkOBCarbonFleeceThick;
   Double_t kLay2       = fgkOBGraphiteFoilThick;
 
-  Double_t xlen, ylen;
   Double_t ymod, zmod;
   Double_t xtru[12], ytru[12];
   Double_t xpos, ypos, ypos1, zpos/*, zpos5cm*/;
-  Double_t zlen;
+  Double_t xlen, ylen, zlen;
   char volname[30];
 
+  Double_t rCoolMin, rCoolMax;
+  if (fStaveModel == AliITSUv1::kOBModel1)
+    rCoolMin = fgkOBCoolTubeInnerDM1/2;
+  else
+    rCoolMin = fgkOBCoolTubeInnerD/2;
+
+  rCoolMax = rCoolMin + fgkOBCoolTubeThick;
 
   zlen = (fNModules*fgkOBModuleZLength + (fNModules-1)*fgkOBModuleGap)/2;
 
@@ -2083,6 +2097,8 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
   TGeoBBox *busKap = new TGeoBBox("BusKap", xHalfSt, fgkOBBusCableKapThick/2,
                                            zlen);
 
+  TGeoBBox *glue = new TGeoBBox("Glue", xHalfSt, fgkOBGlueThick/2, zlen);
+
   TGeoBBox *coldPlate  = new TGeoBBox("ColdPlate", fgkOBHalfStaveWidth/2,
                                      fgkOBColdPlateThick/2, zlen);
 
@@ -2119,11 +2135,15 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
   TGeoBBox *flex2_5cm  = new TGeoBBox("Flex2MV_5cm",xHalfSt,yFlex2/2,flexOverlap/2);
 
   // The half stave container (an XTru to avoid overlaps between neightbours)
+  ylen = ymod + busAl->GetDY() + busKap->GetDY() + coldPlate->GetDY()
+             + graphlat->GetDY() + fleeclat->GetDY();
+  if (fStaveModel == AliITSUv1::kOBModel2)
+    ylen += 2*glue->GetDY();
+
   xtru[0] = xHalfSt;
   ytru[0] = 0;
   xtru[1] = xtru[0];
-  ytru[1] = -2*(ymod + busAl->GetDY() + busKap->GetDY() + coldPlate->GetDY()
-               + graphlat->GetDY() + fleeclat->GetDY());
+  ytru[1] = -2*ylen;
   xtru[2] = fgkOBCoolTubeXDist/2 + fleectub->GetRmax();
   ytru[2] = ytru[1];
   xtru[3] = xtru[2];
@@ -2151,14 +2171,15 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
   TGeoMedium *medCarbonFleece = mgr->GetMedium("ITS_CarbonFleece$"); 
   TGeoMedium *medFGS003       = mgr->GetMedium("ITS_FGS003$"); //amec thermasol
   TGeoMedium *medAir          = mgr->GetMedium("ITS_AIR$");
+  TGeoMedium *medGlue         = mgr->GetMedium("ITS_GLUE$");
 
 
-  TGeoVolume *busAlVol = new TGeoVolume("BusAlVol", busAl , medAluminum);
+  TGeoVolume *busAlVol = new TGeoVolume("PowerBusAlVol", busAl , medAluminum);
   busAlVol->SetLineColor(kCyan);
   busAlVol->SetFillColor(busAlVol->GetLineColor());
   busAlVol->SetFillStyle(4000); // 0% transparent
 
-  TGeoVolume *busKapVol = new TGeoVolume("BusKapVol", busKap, medKapton);
+  TGeoVolume *busKapVol = new TGeoVolume("PowerBusKapVol", busKap, medKapton);
   busKapVol->SetLineColor(kBlue);
   busKapVol->SetFillColor(busKapVol->GetLineColor());
   busKapVol->SetFillStyle(4000); // 0% transparent
@@ -2169,6 +2190,11 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
   coldPlateVol->SetFillColor(coldPlateVol->GetLineColor());
   coldPlateVol->SetFillStyle(4000); // 0% transparent
 
+  TGeoVolume *glueVol = new TGeoVolume("GlueVol", glue, medGlue);
+  glueVol->SetLineColor(kBlack);
+  glueVol->SetFillColor(glueVol->GetLineColor());
+  glueVol->SetFillStyle(4000); // 0% transparent
+
   TGeoVolume *coolTubeVol  = new TGeoVolume("CoolingTubeVol",
                                            coolTube, medKapton);
   coolTubeVol->SetLineColor(kGray);
@@ -2246,36 +2272,61 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
 
   // Now build up the half stave
   ypos = - busKap->GetDY();
-  if (fBuildLevel < 4) // 4 = Kapton
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 4) || // Kapton
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 5) )
     halfStaveVol->AddNode(busKapVol, 1, new TGeoTranslation(0, ypos, 0));
 
   ypos -= (busKap->GetDY() + busAl->GetDY());
-  if (fBuildLevel < 1) // 1 = Aluminum
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 1) || // Aluminum
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 2) )
     halfStaveVol->AddNode(busAlVol, 1, new TGeoTranslation(0, ypos, 0));
 
-  ypos -= (busAl->GetDY() + ymod);
+  ypos -= busAl->GetDY();
+
+  if (fStaveModel == AliITSUv1::kOBModel2) {
+    ypos -= glue->GetDY();
+    if (fBuildLevel < 3) // Glue
+      halfStaveVol->AddNode(glueVol, 1, new TGeoTranslation(0, ypos, 0));
+    ypos -= glue->GetDY();
+  }
+
+  ypos -= ymod;
   for (Int_t j=0; j<fNModules; j++) {
     zpos = -zlen + j*(2*zmod + fgkOBModuleGap) + zmod;
     halfStaveVol->AddNode(moduleVol, j, new TGeoTranslation(0, ypos, zpos));
     fHierarchy[kModule]++;
   }
 
-  ypos -= (ymod + coldPlate->GetDY());
-  if (fBuildLevel < 5) // 5 = Carbon
+  ypos -= ymod;
+
+  if (fStaveModel == AliITSUv1::kOBModel2) {
+    ypos -= glue->GetDY();
+    if (fBuildLevel < 3) // Glue
+      halfStaveVol->AddNode(glueVol, 2, new TGeoTranslation(0, ypos, 0));
+    ypos -= glue->GetDY();
+  }
+
+  ypos -= coldPlate->GetDY();
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 5) || // Carbon
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 6) )
     halfStaveVol->AddNode(coldPlateVol, 1, new TGeoTranslation(0, ypos, 0));
 
   xpos = fgkOBCoolTubeXDist/2;
   ypos1 = ypos - (coldPlate->GetDY() + coolTube->GetRmax());
-  if (fBuildLevel < 3) { // 3 = Water
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 3) || // Water
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 4) ) {
     halfStaveVol->AddNode(coolWaterVol, 1, new TGeoTranslation(-xpos, ypos1, 0));
     halfStaveVol->AddNode(coolWaterVol, 2, new TGeoTranslation( xpos, ypos1, 0));
   }
-  if (fBuildLevel < 4) { // 4 = Kapton
+
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 4) || // Kapton
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 5) ) {
     halfStaveVol->AddNode(coolTubeVol,1, new TGeoTranslation(-xpos, ypos1, 0));
     halfStaveVol->AddNode(coolTubeVol,2, new TGeoTranslation( xpos, ypos1, 0));
   }
 
-  if (fBuildLevel < 5) { // 5 = Carbon
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 5) || // Carbon
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 6) ) {
     halfStaveVol->AddNode(graphtubVol,1, new TGeoTranslation(-xpos, ypos1, 0));
     halfStaveVol->AddNode(graphtubVol,2, new TGeoTranslation( xpos, ypos1, 0));
 
@@ -2285,7 +2336,8 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
 
   xpos = xHalfSt - graphlat->GetDX();
   ypos1 = ypos - (coldPlate->GetDY() + graphlat->GetDY());
-  if (fBuildLevel < 5) { // 5 = Carbon
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 5) || // Carbon
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 6) ) {
     halfStaveVol->AddNode(graphlatVol,1, new TGeoTranslation(-xpos, ypos1, 0));
     halfStaveVol->AddNode(graphlatVol,2, new TGeoTranslation( xpos, ypos1, 0));
 
@@ -2302,7 +2354,8 @@ TGeoVolume* AliITSUv1Layer::CreateStaveModelOuterB1(const TGeoManager *mgr){
 
   xpos = xHalfSt - fleeclat->GetDX();
   ypos1 = ypos - (coldPlate->GetDY() +2*graphlat->GetDY() +fleeclat->GetDY());
-  if (fBuildLevel < 5) { // 5 = Carbon
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 5) || // Carbon
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 6) ) {
     halfStaveVol->AddNode(fleeclatVol,1, new TGeoTranslation(-xpos, ypos1, 0));
     halfStaveVol->AddNode(fleeclatVol,2, new TGeoTranslation( xpos, ypos1, 0));
 
@@ -2371,6 +2424,7 @@ TGeoVolume* AliITSUv1Layer::CreateSpaceFrameOuterB(const TGeoManager *mgr){
       mechStavVol = CreateSpaceFrameOuterBDummy(mgr);
       break;
     case AliITSUv1::kOBModel1:
+    case AliITSUv1::kOBModel2:
       mechStavVol = CreateSpaceFrameOuterB1(mgr);
       break;
     default:
@@ -2703,6 +2757,8 @@ TGeoVolume* AliITSUv1Layer::CreateModuleOuterB(const TGeoManager *mgr){
 //
 // Created:      18 Dec 2013  M. Sitta, A. Barbano
 // Updated:      26 Feb 2014  M. Sitta
+// Updated:      12 Nov 2014  M. Sitta  Model2 is w/o Carbon Plate and Glue
+//                                      and Cu instead of Al
 //
 
 
@@ -2719,7 +2775,7 @@ TGeoVolume* AliITSUv1Layer::CreateModuleOuterB(const TGeoManager *mgr){
 
   // The chip (the same as for IB)
   xlen = (fgkOBHalfStaveWidth/2-xGap/2)/fgkOBNChipRows;
-  ylen = 0.5*fStaveThick; // TO BE CHECKED
+  ylen = 0.5*fChipThick;
   zlen = (fgkOBModuleZLength - (fgkOBChipsPerRow-1)*zGap)/(2*fgkOBChipsPerRow);
 
   TGeoVolume *chipVol = CreateChipInnerB(xlen, ylen, zlen);
@@ -2735,12 +2791,18 @@ TGeoVolume* AliITSUv1Layer::CreateModuleOuterB(const TGeoManager *mgr){
   TGeoBBox *modPlate = new TGeoBBox("CarbonPlate", xlen, ylen, zlen);
 
   // The glue
-  ylen = fgkOBGlueThick/2;
+  ylen = fgkOBGlueThickM1/2;
   TGeoBBox *glue = new TGeoBBox("Glue", xlen, ylen, zlen);
 
-  // The flex cables
-  ylen = fgkOBFlexCableAlThick/2;
-  TGeoBBox *flexAl = new TGeoBBox("FlexAl", xlen, ylen, zlen);
+  // The FPC cables
+  TGeoBBox *flexMetal;
+  if (fStaveModel == AliITSUv1::kOBModel1) {
+    ylen = fgkOBFlexCableAlThick/2;
+    flexMetal = new TGeoBBox("FlexAl", xlen, ylen, zlen);
+  } else {
+    ylen = fgkOBFlexCableCuThick/2;
+    flexMetal = new TGeoBBox("FlexCu", xlen, ylen, zlen);
+  }
 
   ylen = fgkOBFlexCableKapThick/2;
   TGeoBBox *flexKap = new TGeoBBox("FlexKap", xlen, ylen, zlen);
@@ -2748,7 +2810,7 @@ TGeoVolume* AliITSUv1Layer::CreateModuleOuterB(const TGeoManager *mgr){
   // The module
   xlen = fgkOBHalfStaveWidth/2;
   ylen = ychip + modPlate->GetDY() + glue->GetDY() +
-         flexAl->GetDY() + flexKap->GetDY();
+         flexMetal->GetDY() + flexKap->GetDY();
   zlen = fgkOBModuleZLength/2;
   TGeoBBox *module = new TGeoBBox("OBModule", xlen,  ylen, zlen);
 
@@ -2759,6 +2821,7 @@ TGeoVolume* AliITSUv1Layer::CreateModuleOuterB(const TGeoManager *mgr){
   TGeoMedium *medCarbon   = mgr->GetMedium("ITS_CARBON$");
   TGeoMedium *medGlue     = mgr->GetMedium("ITS_GLUE$");
   TGeoMedium *medAluminum = mgr->GetMedium("ITS_ALUMINUM$");
+  TGeoMedium *medCopper   = mgr->GetMedium("ITS_COPPER$");
   TGeoMedium *medKapton   = mgr->GetMedium("ITS_KAPTON(POLYCH2)$");
 
   TGeoVolume *modPlateVol = new TGeoVolume("CarbonPlateVol",
@@ -2772,12 +2835,16 @@ TGeoVolume* AliITSUv1Layer::CreateModuleOuterB(const TGeoManager *mgr){
   glueVol->SetFillColor(glueVol->GetLineColor());
   glueVol->SetFillStyle(4000); // 0% transparent
 
-  TGeoVolume *flexAlVol = new TGeoVolume("FlexAlVol", flexAl, medAluminum);
-  flexAlVol->SetLineColor(kRed);
-  flexAlVol->SetFillColor(flexAlVol->GetLineColor());
-  flexAlVol->SetFillStyle(4000); // 0% transparent
+  TGeoVolume *flexMetalVol;
+  if (fStaveModel == AliITSUv1::kOBModel1)
+    flexMetalVol = new TGeoVolume("FPCAlVol", flexMetal, medAluminum);
+  else
+    flexMetalVol = new TGeoVolume("FPCCuVol", flexMetal, medCopper);
+  flexMetalVol->SetLineColor(kRed);
+  flexMetalVol->SetFillColor(flexMetalVol->GetLineColor());
+  flexMetalVol->SetFillStyle(4000); // 0% transparent
 
-  TGeoVolume *flexKapVol = new TGeoVolume("FlexKapVol", flexKap, medKapton);
+  TGeoVolume *flexKapVol = new TGeoVolume("FPCKapVol", flexKap, medKapton);
   flexKapVol->SetLineColor(kGreen);
   flexKapVol->SetFillColor(flexKapVol->GetLineColor());
   flexKapVol->SetFillStyle(4000); // 0% transparent
@@ -2788,16 +2855,24 @@ TGeoVolume* AliITSUv1Layer::CreateModuleOuterB(const TGeoManager *mgr){
   
 
   // Now build up the module
-  ypos = -module->GetDY() + modPlate->GetDY();
-  if (fBuildLevel < 5) // 5 = Carbon
-    modVol->AddNode(modPlateVol, 1, new TGeoTranslation(0, ypos, 0));
+  // Model1 : CarbonPlate-Glue-Chips-AluminumFPC-KaptonFPC
+  // Model2 : Chips-CopperFPC-KaptonFPCChips
+  ypos = -module->GetDY();
+
+  if (fStaveModel == AliITSUv1::kOBModel1) {
+    ypos += modPlate->GetDY();
+    if (fBuildLevel < 5) // Carbon
+      modVol->AddNode(modPlateVol, 1, new TGeoTranslation(0, ypos, 0));
 
-  ypos += (modPlate->GetDY() + glue->GetDY());
-  if (fBuildLevel < 2) // 2 = Glue
-    modVol->AddNode(glueVol, 1, new TGeoTranslation(0, ypos, 0));
+    ypos += (modPlate->GetDY() + glue->GetDY());
+    if (fBuildLevel < 2) // Glue
+      modVol->AddNode(glueVol, 1, new TGeoTranslation(0, ypos, 0));
+
+    ypos += glue->GetDY();
+  }
 
   xpos = -module->GetDX() + xchip;
-  ypos += (glue->GetDY() + ychip);
+  ypos += ychip;
   for(Int_t k=0; k<fgkOBChipsPerRow; k++)   //put 7x2 chip into one module
     {
       zpos = -module->GetDZ() + zchip + k*(2*zchip + zGap);
@@ -2807,56 +2882,20 @@ TGeoVolume* AliITSUv1Layer::CreateModuleOuterB(const TGeoManager *mgr){
       fHierarchy[kChip]+=2;
     }
 
-  ypos += (ychip + flexAl->GetDY());
-  if (fBuildLevel < 1) // 1 = Aluminum
-    modVol->AddNode(flexAlVol, 1, new TGeoTranslation(0, ypos, 0));
+  ypos += (ychip + flexMetal->GetDY());
+  if (fBuildLevel < 1) // Model1: Aluminum  Model2: Copper
+    modVol->AddNode(flexMetalVol, 1, new TGeoTranslation(0, ypos, 0));
 
-  ypos += (flexAl->GetDY() + flexKap->GetDY());
-  if (fBuildLevel < 4) // 4 = Kapton
+  ypos += (flexMetal->GetDY() + flexKap->GetDY());
+  if ( (fStaveModel == AliITSUv1::kOBModel1 && fBuildLevel < 4) || // Kapton
+       (fStaveModel == AliITSUv1::kOBModel2 && fBuildLevel < 5) )
     modVol->AddNode(flexKapVol, 1, new TGeoTranslation(0, ypos, 0));
-  
 
   // Done, return the module
   return modVol;
 }
 
 //________________________________________________________________________
-Double_t AliITSUv1Layer::RadiusOfTurboContainer(){
-//
-// Computes the inner radius of the air container for the Turbo configuration
-// as the radius of either the circle tangent to the stave or the circle
-// passing for the stave's lower vertex
-//
-// Input:
-//         none (all needed parameters are class members)
-//
-// Output:
-//
-// Return:
-//        the radius of the container if >0, else flag to use the lower vertex
-//
-// Created:      08 Mar 2012  Mario Sitta
-//
-
-  Double_t rr, delta, z, lstav, rstav;
-
-  if (fStaveThick > 89.) // Very big angle: avoid overflows since surely
-    return -1;            // the radius from lower vertex is the right value
-
-  rstav = fLayRadius + 0.5*fStaveThick;
-  delta = (0.5*fStaveThick)/CosD(fStaveTilt);
-  z     = (0.5*fStaveThick)*TanD(fStaveTilt);
-
-  rr = rstav - delta;
-  lstav = (0.5*fStaveWidth) - z;
-
-  if ( (rr*SinD(fStaveTilt) < lstav) )
-    return (rr*CosD(fStaveTilt));
-  else
-    return -1;
-}
-
-//________________________________________________________________________
 void AliITSUv1Layer::SetNUnits(Int_t u)
 {
 //
index a021f11..9491a8c 100644 (file)
@@ -40,11 +40,11 @@ class AliITSUv1Layer : public AliITSv11Geometry {
     //
     Bool_t    IsTurbo() const {return fIsTurbo;};
 
-    Double_t  GetStaveThick() const {return fStaveThick;};
-    Double_t  GetStaveTilt()  const {return fStaveTilt;};
-    Double_t  GetStaveWidth() const {return fStaveWidth;};
+    Double_t  GetChipThick()   const {return fChipThick;};
+    Double_t  GetStaveTilt()   const {return fStaveTilt;};
+    Double_t  GetStaveWidth()  const {return fStaveWidth;};
     Double_t  GetSensorThick() const {return fSensorThick;};
-    Double_t  GetNStaves()    const {return fNStaves;};
+    Double_t  GetNStaves()     const {return fNStaves;};
     Double_t  GetNChips()      const {return fNChips;};
     Double_t  GetRadius()      const {return fLayRadius;};
     Double_t  GetPhi0()        const {return fPhi0;};
@@ -56,9 +56,10 @@ class AliITSUv1Layer : public AliITSv11Geometry {
     Int_t     GetNModulesPerParent()    const {return fHierarchy[kModule];}
     Int_t     GetNChipsPerParent()      const {return fHierarchy[kChip];}
     //
+    Int_t     GetBuildLevel()           const {return fBuildLevel;}
     AliITSUv1::AliITSUModel_t GetStaveModel() const {return fStaveModel;}
     //
-    void      SetStaveThick(Double_t t)      {fStaveThick = t;};
+    void      SetChipThick(Double_t t)       {fChipThick = t;};
     void      SetStaveTilt(Double_t t);
     void      SetStaveWidth(Double_t w);
     void      SetSensorThick(Double_t t)     {fSensorThick = t;};
@@ -75,8 +76,6 @@ class AliITSUv1Layer : public AliITSv11Geometry {
   private:
     void CreateLayerTurbo(TGeoVolume *moth);
 
-    Double_t RadiusOfTurboContainer();
-
     TGeoVolume* CreateStave(const TGeoManager *mgr=gGeoManager);
     //TGeoVolume* CreateChip(Double_t x, Double_t z, const TGeoManager *mgr=gGeoManager);
     TGeoVolume* CreateModuleInnerB(Double_t x,Double_t y, Double_t z, const TGeoManager *mgr=gGeoManager);
@@ -96,7 +95,7 @@ class AliITSUv1Layer : public AliITSv11Geometry {
     TGeoVolume* CreateStaveOuterB(const TGeoManager *mgr=gGeoManager);
     TGeoVolume* CreateStaveModelOuterBDummy(const TGeoManager *mgr=gGeoManager) const;
     TGeoVolume* CreateStaveModelOuterB0(const TGeoManager *mgr=gGeoManager);
-    TGeoVolume* CreateStaveModelOuterB1(const TGeoManager *mgr=gGeoManager);
+    TGeoVolume* CreateStaveModelOuterB12(const TGeoManager *mgr=gGeoManager);
     TGeoVolume* CreateSpaceFrameOuterB(const TGeoManager *mgr=gGeoManager);
     TGeoVolume* CreateSpaceFrameOuterBDummy(const TGeoManager *mgr=gGeoManager) const;
     TGeoVolume* CreateSpaceFrameOuterB1(const TGeoManager *mgr=gGeoManager);
@@ -117,7 +116,7 @@ class AliITSUv1Layer : public AliITSv11Geometry {
     Double_t  fLayRadius;   // Inner radius of this layer
     Double_t  fZLength;     // Z length of this layer
     Double_t  fSensorThick; // Sensor thickness
-    Double_t  fStaveThick;  // Stave thickness
+    Double_t  fChipThick;   // Chip thickness
     Double_t  fStaveWidth;  // Stave width (for turbo layers only)
     Double_t  fStaveTilt;   // Stave tilt angle (for turbo layers only) in degrees
     Int_t     fNStaves;     // Number of staves in this layer
@@ -137,7 +136,7 @@ class AliITSUv1Layer : public AliITSv11Geometry {
     static const Int_t    fgkNumberOfInnerLayers;// Number of IB Layers
 
     static const Double_t fgkDefaultSensorThick; // Default sensor thickness
-    static const Double_t fgkDefaultStaveThick;  // Default stave thickness
+    static const Double_t fgkDefaultChipThick;   // Default chip thickness
 
     // Inner Barrel Parameters
     static const Int_t    fgkIBChipsPerRow;      // IB chips per row in module
@@ -153,16 +152,20 @@ class AliITSUv1Layer : public AliITSv11Geometry {
     static const Double_t fgkOBChipXGap;         // Gap between OB chips on X
     static const Double_t fgkOBChipZGap;         // Gap between OB chips on Z
     static const Double_t fgkOBFlexCableAlThick; // Thickness of FPC Aluminum
+    static const Double_t fgkOBFlexCableCuThick; // Thickness of FPC Copper
     static const Double_t fgkOBFlexCableKapThick;// Thickness of FPC Kapton
     static const Double_t fgkOBBusCableAlThick;  // Thickness of Bus Aluminum
     static const Double_t fgkOBBusCableKapThick; // Thickness of Bus Kapton
     static const Double_t fgkOBCarbonPlateThick; // OB Carbon Plate Thickness
     static const Double_t fgkOBColdPlateThick;   // OB Cold Plate Thickness
-    static const Double_t fgkOBGlueThick;        // OB Glue total Thickness
+    static const Double_t fgkOBGlueThickM1;      // OB Glue total Thickness
+    static const Double_t fgkOBGlueThick;        // OB Glue Thickness in Model2
     static const Double_t fgkOBModuleZLength;    // OB Chip Length along Z
     static const Double_t fgkOBHalfStaveYTrans;  // OB half staves Y transl.
     static const Double_t fgkOBHalfStaveXOverlap;// OB half staves X overlap
     static const Double_t fgkOBGraphiteFoilThick;// OB graphite foil thickness
+    static const Double_t fgkOBCarbonFleeceThick;// OB carbon fleece thickness
+    static const Double_t fgkOBCoolTubeInnerDM1; // OB cooling inner diameter
     static const Double_t fgkOBCoolTubeInnerD;   // OB cooling inner diameter
     static const Double_t fgkOBCoolTubeThick;    // OB cooling tube thickness
     static const Double_t fgkOBCoolTubeXDist;    // OB cooling tube separation
index 4aee313..9dfe298 100644 (file)
@@ -31,8 +31,8 @@ void CreateITSUv1()
   const double kPitchZ = 20e-4;
   const int    kNRow   = 650; 
   const int    kNCol   = 1500;
-  const double kSiThickIB = 150e-4;
-  const double kSiThickOB = 150e-4;
+  const double kSiThickIB = 50e-4;
+  const double kSiThickOB = 50e-4;
   //  const double kSensThick = 120e-4;   // -> sensor Si thickness
   //
   const double kReadOutEdge = 0.2;   // width of the readout edge (passive bottom)
@@ -78,7 +78,7 @@ void CreateITSUv1()
   //
   AliITSUv1 *ITS  = new AliITSUv1("ITS Upgrade",kNLr);
   ITS->SetStaveModelIB(AliITSUv1::kIBModel22);
-  ITS->SetStaveModelOB(AliITSUv1::kOBModel1);
+  ITS->SetStaveModelOB(AliITSUv1::kOBModel2);
   //
   const int kNWrapVol = 3;
   const double wrpRMin[kNWrapVol]  = { 2.1, 15.0, 32.0};
diff --git a/ITS/UPGRADE/testITSUv1/CreateITSUv1.C_template b/ITS/UPGRADE/testITSUv1/CreateITSUv1.C_template
new file mode 100644 (file)
index 0000000..b2855c6
--- /dev/null
@@ -0,0 +1,116 @@
+#if !defined(__CINT__) || defined(__MAKECINT__)
+#include <TSystem.h>
+#include <TMath.h>
+#endif
+
+//---------------------------------------
+double radii2Turbo(double rMin,double rMid,double rMax, double sensW)
+{
+  // compute turbo angle from radii and sensor width
+  return TMath::ASin((rMax*rMax-rMin*rMin)/(2*rMid*sensW))*TMath::RadToDeg();
+}
+
+double radii2Phi(double rMin,double rMid,double rMax, double sensW)
+{
+  // compute phi coverage
+  return 2*TMath::ACos((rMax+rMin)*
+                      (rMid*rMid+rMin*rMax-sensW*sensW/4.)/
+                      (4.*rMid*rMax*rMin));
+}
+
+void CreateITSUv1()
+{
+  //
+  gSystem->Load("libITSUpgradeBase.so");
+  gSystem->Load("libITSUpgradeSim.so");
+  //
+  // build ITS upgrade detector
+  // sensitive area 13x15mm (X,Z) with 20x20 micron pitch, 2mm dead zone on readout side and 50 micron guardring
+  const double kSensThick = 18e-4;
+  const double kPitchX = 20e-4;
+  const double kPitchZ = 20e-4;
+  const int    kNRow   = 650; 
+  const int    kNCol   = 1500;
+//  const double kSiThickIB = 250e-4;
+  const double kSiThickIB = 50e-4;
+  const double kSiThickOB = 50e-4;
+  //  const double kSensThick = 120e-4;   // -> sensor Si thickness
+  //
+  const double kReadOutEdge = 0.2;   // width of the readout edge (passive bottom)
+  const double kGuardRing   = 50e-4; // width of passive area on left/right/top of the sensor
+  //
+  const int kNLr = 7;
+  const int kNLrInner = 3;
+  const int kBuildLevel = BUILDLEVEL;
+  enum {kRmn,kRmd,kRmx,kNModPerStave,kPhi0,kNStave,kNPar};
+  // Radii are from last TDR (ALICE-TDR-017.pdf Tab. 1.1, rMid is mean value)
+  const double tdr5dat[kNLr][kNPar] = { 
+    {2.24, 2.34, 2.67,  9., 16.37, 12}, // for each inner layer: rMin,rMid,rMax,NChip/Stave, phi0, nStaves
+    {3.01, 3.15, 3.46,  9., 12.03, 16},
+    {3.78, 3.93, 4.21,  9., 10.02, 20},
+    {-1,  19.6 ,   -1,  4.,  0.  , 24},  // for others: -, rMid, -, NMod/HStave, phi0, nStaves // 24 was 49
+    {-1,  24.55, -1,    4.,  0.  , 30},  // 30 was 61
+    {-1,  34.39, -1,    7.,  0.  , 42},  // 42 was 88
+    {-1,  39.34, -1,    7.,  0.  , 48}   // 48 was 100
+  };
+  const int nChipsPerModule = 7; // For OB: how many chips in a row
+
+  // create segmentations:
+  AliITSUSegmentationPix* seg0 = new AliITSUSegmentationPix(0,        // segID (0:9)
+                                                           1,  // chips per module
+                                                           kNCol,    // ncols (total for module)
+                                                           kNRow,    // nrows
+                                                           kPitchX,  // default row pitch in cm
+                                                           kPitchZ,  // default col pitch in cm
+                                                           kSensThick,  // sensor thickness in cm
+                                                           -1,     // no special left col between chips
+                                                           -1,     // no special right col between chips
+                                                           kGuardRing, // left
+                                                           kGuardRing, // right
+                                                           kGuardRing, // top
+                                                           kReadOutEdge  // bottom
+                                                           );    // see AliITSUSegmentationPix.h for extra options
+  seg0->Store(AliITSUGeomTGeo::GetITSsegmentationFileName());
+  //
+  seg0->Print();
+  //
+  double dzLr,rLr,phi0,turbo;
+  int nStaveLr,nModPerStaveLr,idLr;
+  //
+  AliITSUv1 *ITS  = new AliITSUv1("ITS Upgrade",kNLr);
+  ITS->SetStaveModelIB(AliITSUv1::kIBModel22);
+  ITS->SetStaveModelOB(AliITSUv1::kOBModel2);
+  //
+  const int kNWrapVol = 3;
+  const double wrpRMin[kNWrapVol]  = { 2.1, 15.0, 32.0};
+  const double wrpRMax[kNWrapVol]  = { 7.0, 27.0+2.5, 43.0+1.5};
+  const double wrpZSpan[kNWrapVol] = {28.0, 86.0, 150.0};
+  //
+  ITS->SetNWrapVolumes(kNWrapVol); // define wrapper volumes for layers
+  for (int iw=0;iw<kNWrapVol;iw++) ITS->DefineWrapVolume(iw,wrpRMin[iw],wrpRMax[iw],wrpZSpan[iw]);
+  //
+  for (int idLr=0;idLr<kNLr;idLr++) {
+    rLr   = tdr5dat[idLr][kRmd];
+    phi0  = tdr5dat[idLr][kPhi0]; 
+    //
+    nStaveLr = TMath::Nint(tdr5dat[idLr][kNStave]);
+    nModPerStaveLr =  TMath::Nint(tdr5dat[idLr][kNModPerStave]);
+    int nChipsPerStaveLr = nModPerStaveLr;
+    //
+    if (idLr>=kNLrInner) {
+      nChipsPerStaveLr *= nChipsPerModule;
+      ITS->DefineLayer(idLr, phi0, rLr, nChipsPerStaveLr*seg0->Dz(), nStaveLr, nModPerStaveLr, 
+                      kSiThickOB, seg0->Dy(), seg0->GetChipTypeID(),kBuildLevel);
+      //      printf("Add Lr%d: R=%6.2f DZ:%6.2f Staves:%3d NMod/Stave:%3d\n",
+      //            idLr,rLr,nChipsPerStaveLr*seg0->Dz(),nStaveLr,nModPerStaveLr);
+    } else {
+      turbo = -radii2Turbo(tdr5dat[idLr][kRmn],rLr,tdr5dat[idLr][kRmx],seg0->Dx());    
+      ITS->DefineLayerTurbo(idLr, phi0, rLr, nChipsPerStaveLr*seg0->Dz(), nStaveLr, nChipsPerStaveLr, 
+                           seg0->Dx(), turbo, kSiThickIB, seg0->Dy(), seg0->GetChipTypeID(),kBuildLevel);
+      //      printf("Add Lr%d: R=%6.2f DZ:%6.2f Turbo:%+6.2f Staves:%3d NMod/Stave:%3d\n",
+      //            idLr,rLr,nChipsPerStaveLr*seg0->Dz(),turbo,nStaveLr,nModPerStaveLr);
+    }
+    //
+  }
+  //  
+}
diff --git a/ITS/UPGRADE/testITSUv1/GetMaterialBudget.C b/ITS/UPGRADE/testITSUv1/GetMaterialBudget.C
new file mode 100644 (file)
index 0000000..fe6a777
--- /dev/null
@@ -0,0 +1,760 @@
+#ifndef __CINT__
+
+#include "TSystem.h"
+#include "TGeoManager.h"
+#include <AliRunLoader.h>
+#include <AliGeomManager.h>
+#include <AliITSgeom.h>
+#include <AliTracker.h>
+#include "TRandom.h"
+#include "TMath.h"
+#include "TH2F.h"
+#include "TCanvas.h"
+#include "TFile.h"
+#include <AliRun.h>
+#include <TLegend.h>
+
+#include <TGeoVolume.h>
+#include "TGeoMaterial.h"
+#include "TGeoMedium.h"
+#include <TGeoTube.h>
+#include <AliITSUGeomTGeo.h>
+#include <TPaveText.h>
+#include <TText.h>
+
+#endif
+
+/*
+  gSystem->Load("libITSUpgradeBase");
+  gSystem->Load("libITSUpgradeSim");
+
+  .L GetMaterialBudget.C
+
+  DrawMaterialBudget_SPLITTED(); 
+  DrawMaterialBudget_FromTo();
+
+  GetMaterialBudget_EtaVsPhi(0,0);
+  GetMaterialBudget_inPhi(0,0);
+
+*/
+
+// Original version
+// Chinorat Kobdaj (kobdaj@g.sut.ac.th)
+// Revised and adapted to newer AliITSUv1 versions
+// Mario Sitta <sitta@to.infn.it> - Nov 2014
+
+
+enum {
+    kIBModelDummy=0,
+    kIBModel0=1,
+    kIBModel1=2, 
+    kIBModel21=3,
+    kIBModel22=4,
+    kIBModel3=5,
+    kOBModelDummy=6,
+    kOBModel0=7,
+    kOBModel1=8, 
+    kOBModel2=9 
+};
+
+
+//______________________________________________________________________
+void PrintMaterialDefs(const char *geofile="geometry.root")
+{
+  //
+  // Print material definition for all ITS materials/mediums
+  // Rewritten - M.Sitta 15 Nov 2014
+  //
+  TGeoManager::Import(geofile);
+
+  TGeoMedium *med;
+  TGeoMaterial *mat;
+  char mediaName[50], matName[50], shortName[50];
+
+  Int_t nMedia = gGeoManager->GetListOfMedia()->GetEntries();
+
+  printf("\n\n\n");
+  printf("              ==== ITSU  Material  Properties ====\n\n");
+  printf("    A      Z   d (g/cm3)  RadLen (cm)  IntLen (cm)\t Name\n");
+
+  // Loop on media, select ITS materials, print their characteristics
+  for (Int_t i = 0; i < nMedia; i++) {
+    med = (TGeoMedium *)(gGeoManager->GetListOfMedia()->At(i));
+    strcpy(mediaName, med->GetName());
+    // The trick here: the name must begin with the string "ITS_"
+    // so the pointer to the first occurrence of the substring as returned
+    // by strstr() must match the beginning of the string
+    if (strstr(mediaName,"ITS_") == mediaName) { // name begins with "ITS_"
+      mat = med->GetMaterial();
+      strcpy(matName, mat->GetName());
+      strncpy(shortName, matName+4, strlen(matName)-5); // get rid of ITS_ , $
+      shortName[strlen(matName)-5] = '\0';
+      printf(" %5.1f %6.1f %8.3f %13.3f %11.3f\t %s\n", 
+            mat->GetA(), mat->GetZ(), mat->GetDensity(), mat->GetRadLen(),
+            mat->GetIntLen(), shortName);
+    }
+  }
+
+  printf("\n");
+
+//   printf("Values in Corrados estimate");
+//   const Int_t nC=6;
+//   TString medNC[nC]={"ITS_SI$","ITS_WATER$","ITS_CARBON$","ITS_KAPTON$","ITS_FLEXCABLE$","ITS_GLUE$"};
+//   Double_t radl[nC]={9.36,36.1,25,28.6,13.3,44.37};
+//   for (Int_t i=0; i<nC; i++) {
+//     printf("%8.2lf (cm)  \t%s\n", radl[i],medNC[i].Data());
+//   }
+  
+}
+
+
+//______________________________________________________________________
+Double_t ComputeMaterialBudget(Double_t rmin, Double_t rmax, Double_t etaPos,
+                              Double_t phiMax,
+                              TH1F* xOverX01d, TH2F* xOverX02d = 0)
+{
+  //
+  // Ancillary function to compute material budget between rmin and rmax
+  // and +/- etaPos in the 0--phiMax range, filling two histograms
+  // Rewritten - M.Sitta 15 Nov 2014
+  //
+  Int_t n = 5e4; //testtracks
+
+  TH1F *nParticles1d = new TH1F("", "", 300, 0, phiMax);
+  TH2F *nParticles2d = new TH2F("", "", 300, -etaPos, etaPos, 300, 0, phiMax);
+
+  Double_t mparam[7]={0.,0.,0.,0.,0.};
+  Double_t point1[3],point2[3];
+
+  for (Int_t it=0; it<n; it++) {
+
+    // PHI VS ETA
+    Double_t phi = gRandom->Rndm()*phiMax;
+    Double_t eta = gRandom->Rndm()*2*etaPos - etaPos;
+    Double_t theta = TMath::ATan(TMath::Exp(-eta))*2;
+    point1[0] = rmin*TMath::Cos(phi);
+    point1[1] = rmin*TMath::Sin(phi);
+    point1[2] = rmin/TMath::Tan(theta);
+    point2[0] = rmax*TMath::Cos(phi);
+    point2[1] = rmax*TMath::Sin(phi);
+    point2[2] = rmax/TMath::Tan(theta);
+    
+    AliTracker::MeanMaterialBudget(point1,point2,mparam);
+
+    xOverX01d->Fill(phi, mparam[1]*100);
+    nParticles1d->Fill(phi, 1.);
+
+    if (xOverX02d != 0) {
+      xOverX02d->Fill(eta, phi, mparam[1]*100);
+      nParticles2d->Fill(eta, phi, 1.);
+    }
+
+    if (!(it%10000)) cout<<" : "<<mparam[1]*100<<"   phi:"<<phi<<endl;
+
+  }
+
+  // normalization to number of particles in case of phi vs eta
+  Double_t theta = TMath::ATan(TMath::Exp(-etaPos/2))*2;
+  printf("<eta>=%lf -> Sin(theta) %lf\n",etaPos/2,TMath::Sin(theta)); 
+
+  for (Int_t ix = 1; ix<=nParticles1d->GetNbinsX(); ix++) {
+    if (nParticles1d->GetBinContent(ix) > 0) 
+      xOverX01d->SetBinContent(ix,xOverX01d->GetBinContent(ix)/nParticles1d->GetBinContent(ix)*TMath::Sin(theta));
+    }
+
+  if (xOverX02d) {
+    for (Int_t ix = 1; ix<=nParticles2d->GetNbinsX(); ix++) {
+      for (Int_t iy = 1; iy<=nParticles2d->GetNbinsY(); iy++) {
+       if (nParticles2d->GetBinContent(ix,iy) > 0) 
+         xOverX02d->SetBinContent(ix,iy,xOverX02d->GetBinContent(ix,iy)/nParticles2d->GetBinContent(ix,iy));
+      }
+    }
+  }
+
+  Double_t mean = 0;
+  for (Int_t ix = 1; ix<=xOverX01d->GetNbinsX(); ix++)
+    mean+=xOverX01d->GetBinContent(ix);
+
+  mean /= xOverX01d->GetNbinsX();
+
+  return mean;
+}
+
+
+//______________________________________________________________________
+void DrawMaterialBudget_Splitted(Int_t nLay = 0, Double_t rmin = 1.,
+                                Double_t rmax = 5.)
+{
+  //
+  // Function to factorize the percentage of each medium in the
+  // total material budget between rmin and rmax for layer nLay
+  // Cleaned up - M.Sitta 15 Nov 2014
+  //
+  TCanvas *c2 = new TCanvas("c2","c2");
+
+  Double_t etaPos = 1.;
+
+  Bool_t firstPlot = 1;
+
+  // Check parameters
+  if (nLay < 0 || nLay > 6) {
+    printf("ERROR! Wrong layer number %d\n",nLay);
+    return;
+  }
+
+  if (rmin < 0. || rmax < 0. || rmax < rmin) {
+    printf("ERROR! Wrong radial values rmin = %f rmax = %f\n",rmin,rmax);
+    return;
+  }
+
+  // In current version, build levels have different meaning for IB and OB
+  const Int_t nScenariosIB = 6, nScenariosOB = 7;
+  TString strDescripOB[nScenariosOB] = {"Copper", "Aluminum", "Glue", "Water",
+                                       "Kapton", "Carbon"  , "Silicon"      };
+  //  Color codes
+  //  Water                kBlack
+  //  Kapton               kYellow
+  //  Carbon               kRed
+  //  Glue                 kBlue
+  //  Aluminum             kCyan
+  //  Si-Sensor            kGreen
+  //  Copper               kGray
+  Int_t colorsOB[nScenariosOB] = {kGray  , kCyan, kBlue, kBlack,
+                                 kYellow, kRed , kGreen       };
+
+//   TString strDescripIB[nScenariosIB] = {"Carbon Structure", "Water",
+//                                     "Cooling Pipe Walls and ColdPlate",
+//                                     "Glue", "Flex Cable", "Pixel Chip"};
+  TString strDescripIB[nScenariosIB] = {"Flex cable", "Glue",
+                                       "Carbon structure", "Water",
+                                       "Cooling walls", "Pixel Chip"};
+  //  Color codes
+  //  Flex Cable           kBlack
+  //  Glue                 kYellow
+  //  Carbon Structure     kRed
+  //  Water                kBlue
+  //  Cooling Walls        kCyan 
+  //  Si-Sensor            kGreen 
+  Int_t colorsIB[nScenariosIB] = {kCyan, kBlue, kBlack, kYellow, kRed, kGreen};
+
+  // Setup
+  const Int_t maxNumberOfScenarios = nScenariosOB;
+
+  Double_t meanX0[maxNumberOfScenarios];
+  TH1F *l[maxNumberOfScenarios+1];
+
+  TString strDescrip[maxNumberOfScenarios];
+  Int_t colors[maxNumberOfScenarios];
+
+  // Choose which scenario based on Layer number and Model
+  Int_t nScenarios;
+
+  TGeoManager::Import(Form("geometry_0.root"));
+  AliITSUGeomTGeo* gm = new AliITSUGeomTGeo(kTRUE);
+  Int_t nLad = gm->GetNStaves(nLay);
+
+  Float_t phiMax = TMath::TwoPi()/nLad*2;
+
+  char title[30];
+  Int_t model, buildlevel;
+  if (nLay < 3) { // IB has only 6 scenarios
+    nScenarios = nScenariosIB;
+    for (Int_t i=0; i<nScenarios; i++) {
+      strDescrip[i] = strDescripIB[i];
+      colors[i] = colorsIB[i];
+    }
+  } else {
+    strcpy(title,
+          gGeoManager->GetVolume(Form("ITSULayer%d",nLay))->GetTitle());
+    if (strlen(title) == 0) { // Old OB has only 6 scenarios like IB
+      nScenarios = nScenariosIB;
+      for (Int_t i=0; i<nScenarios; i++) {
+       strDescrip[i] = strDescripOB[i+1];
+       colors[i] = colorsOB[i+1];
+      }
+    } else { // New OB: check model
+      sscanf(title, "Model %d Build level %d", &model, &buildlevel);
+      if (model == kOBModel2) {
+       nScenarios = nScenariosOB;
+       for (Int_t i=0; i<nScenarios; i++) {
+         strDescrip[i] = strDescripOB[i];
+         colors[i] = colorsOB[i];
+       }
+      } else {
+       nScenarios = nScenariosIB;
+       for (Int_t i=0; i<nScenarios; i++) {
+         strDescrip[i] = strDescripOB[i+1];
+         colors[i] = colorsOB[i+1];
+       }
+      }
+    }
+  } // if (nLay < 3)
+
+  delete gGeoManager;
+
+  for (Int_t i=0; i<nScenarios; i++) {
+
+    printf(" -> Loading geometry_%d.root .....\n",i);
+    TGeoManager::Import(Form("geometry_%d.root",i)); 
+
+    strcpy(title,
+          gGeoManager->GetVolume(Form("ITSULayer%d",nLay))->GetTitle());
+    if (strlen(title) != 0) {
+      sscanf(title, "Model %d Build level %d", &model, &buildlevel);
+      if (i != buildlevel)
+       printf("WARNING! Possible mismatch: file geometry_%d.root created with Build level %d\n",i,buildlevel);
+    }
+
+    TH1F *xOverX01d =  new TH1F("", "", 300, 0, phiMax);
+
+    Double_t mean = ComputeMaterialBudget(rmin, rmax, etaPos, phiMax,
+                                         xOverX01d);
+    meanX0[i] = mean;
+    cout<<"Mean X/X0: " << meanX0[i] << " %" << endl;
+
+    xOverX01d->GetXaxis()->SetTitle("#phi (rad)");
+    xOverX01d->GetYaxis()->SetTitle(Form("X/X_{0} (%%) at #eta=0 "));
+    xOverX01d->SetFillColor(colors[i]-7);
+    if (i==0)  xOverX01d->SetFillColor(colors[i]);
+    if ( (nScenarios == nScenariosIB && i==2) ||
+        (nScenarios == nScenariosOB && i==3) )
+      xOverX01d->SetFillColor(colors[i]);
+    xOverX01d->SetStats(0);
+
+    l[i] = new TH1F("","",1,0,phiMax);
+    l[i]->SetBinContent(1,mean);
+    l[i]->SetStats(0);
+    l[i]->SetLineColor(colors[i]-7);
+    if (i==0) l[i]->SetLineColor(kGray);
+    if ( (nScenarios == nScenariosIB && i==2) ||
+        (nScenarios == nScenariosOB && i==3) ) l[i]->SetLineColor(12);
+
+    c2->cd();
+    if (firstPlot) {
+      xOverX01d->SetMinimum(0.);
+      xOverX01d->DrawCopy(); 
+      firstPlot=0;
+    } else
+      xOverX01d->DrawCopy("same"); 
+     
+    delete gGeoManager;
+  }
+  
+  // Build meaningful legend
+  TLegend *leg = new TLegend(0.76,0.77,0.99,0.99,"");
+  leg->SetFillColor(0);
+
+  for (Int_t i=0; i<nScenarios; i++) {
+    // contribution in percent
+    Double_t contr = 0;
+    if (i == nScenarios-1)
+      contr = (meanX0[i])/meanX0[0]*100; 
+    else
+      contr = (meanX0[i]-meanX0[i+1])/meanX0[0]*100; 
+    strDescrip[i].Append(Form(" (%3.1lf%%)",contr));
+    leg->AddEntry(l[i],strDescrip[i].Data(),"l");
+  }
+
+  TPaveText *pt = new TPaveText(0.76,0.70,0.99,0.77,"brNDC");
+  pt->SetBorderSize(1); // no shadow
+  pt->SetTextFont(12);
+  pt->SetFillColor(0);
+  pt->AddText(Form("Mean X/X0 = %4.3lf%%",meanX0[0]));
+  pt->Draw();
+
+  leg->Draw();
+
+  l[nScenarios]=(TH1F*)l[0]->Clone();
+  l[nScenarios]->SetLineColor(1);
+  l[nScenarios]->SetLineWidth(2);
+  l[nScenarios]->DrawCopy("same");
+  
+  c2->SaveAs("Material-details.pdf");
+  
+}
+
+
+//______________________________________________________________________
+void DrawMaterialBudget_FromTo(Int_t nLay = 0, Double_t rmin = 1.,
+                              Double_t rmax = 5., Bool_t only2D = 1)
+{
+  //
+  // Function to compute material budget between rmin and rmax
+  // for layer nLay. If only2D is false a 1D histogram is also plotted
+  // Cleaned up and simplified - M.Sitta 18 Nov 2014
+  //
+
+  Double_t etaPos = 0.25;
+
+  // Check parameters
+  if (nLay < 0 || nLay > 6) {
+    printf("ERROR! Wrong layer number %d\n",nLay);
+    return;
+  }
+
+  if (rmin < 0. || rmax < 0. || rmax < rmin) {
+    printf("ERROR! Wrong radial values rmin = %f rmax = %f\n",rmin,rmax);
+    return;
+  }
+
+
+  TGeoManager::Import("geometry.root"); 
+  AliITSUGeomTGeo* gm = new AliITSUGeomTGeo(kTRUE);
+  Int_t nLad = gm->GetNStaves(nLay);
+
+  Float_t phiMax = TMath::TwoPi()/nLad*2;
+
+  TH2F *xOverX0   =  new TH2F("", "", 300, -1, 1., 300, 0, phiMax);
+  TH1F *xOverX01d =  new TH1F("", "", 300,  0, phiMax);
+
+  Double_t mean = ComputeMaterialBudget(rmin, rmax, -1, phiMax,
+                                       xOverX01d, xOverX0);
+  cout<<"Mean X/X0: " << mean << " %" << endl;
+
+  TH1F *l = new TH1F("", "", 1, 0, phiMax);
+  l->SetBinContent(1, mean);
+  l->SetLineColor(2);
+  l->SetStats(0);
+  l->SetLineWidth(2);
+
+  // Draw the histograms
+  xOverX0->SetTitle(0);
+  xOverX0->GetXaxis()->SetTitle("pseudorapidity #eta ");
+  xOverX0->GetYaxis()->SetTitle("#phi (rad)");
+  xOverX0->GetZaxis()->SetTitle("X/X_{0} (%)");
+  xOverX0->SetStats(0);
+
+  xOverX01d->SetTitle(Form("X/X_{0} (%%) within r=(%2.1lf-%2.1lf)cm average over #eta=(%1.1lf,%1.1lf)", rmin, rmax, -1., 1.));
+  xOverX01d->GetXaxis()->SetTitle("#phi (rad)");
+  xOverX01d->GetYaxis()->SetTitle("X/X_{0} (%) average over #eta=(-1,1)");
+  xOverX01d->SetStats(0);
+
+  if (!only2D) {
+    TCanvas *c1 = new TCanvas("c1","c1",800,800);
+    c1->Divide(1,2);
+    c1->cd(1); xOverX0->Draw("colz");
+    c1->cd(2); xOverX01d->DrawCopy(); l->DrawCopy("same");
+    c1->SaveAs("Material-1D.pdf");
+  } else {
+    TCanvas *c1 = new TCanvas("c1","c1");
+    xOverX0->Draw("colz");
+    c1->SaveAs("Material-2D.pdf");
+  }
+
+}
+
+
+//______________________________________________________________________
+void ComputeGeometricalParameters(const char *geofile, Double_t *rmin,
+                                 Double_t *rmax, Double_t *zPos,
+                                 Int_t *nlad, Bool_t fromGDML = 0)
+{
+  //
+  // Ancillary function to retrieve various geometrical parameters
+  // from a geometry file. The caller must ensure that all vectors
+  // have proper dimensions: since this function is designed to be used
+  // internally, no check is made on parameters!
+  // Rewritten - M.Sitta 20 Nov 2014
+  //
+
+  if (fromGDML)  // from GDML
+    TGeoManager::Import(geofile);
+//    gGeoManager->ViewLeaves(true); 
+//    gGeoManager->GetTopVolume()->Draw("ogl"); 
+  else {         // from AliRoot simulation
+        
+//    gAlice=NULL;
+//    AliRunLoader* runLoader = AliRunLoader::Open("galice.root");
+//    runLoader->LoadgAlice();
+    
+//    gAlice = runLoader->GetAliRun();
+//    AliGeomManager::LoadGeometry(geofile);
+    TGeoManager::Import(geofile);
+    AliITSUGeomTGeo* gm = new AliITSUGeomTGeo(kTRUE);
+
+    TGeoVolume *pipeV = gGeoManager->GetVolume("IP_PIPE");
+    if (!pipeV) {
+      printf("Pipe volume %s is not in the geometry\n", "IP_PIPE");
+      return;
+    } else {
+      printf("%s   ","IP_PIPE");
+      TGeoTube *t = (TGeoTube*)(pipeV->GetShape());
+      printf(" r = (%3.2lf %3.2lf) ", t->GetRmin(), t->GetRmax());
+      printf(" z = %3.2lf\n", t->GetDz());
+
+      rmin[0] = t->GetRmin();
+      rmax[0] = t->GetRmax();
+      nlad[0] = 0;
+      zPos[0] = 0;
+    }
+
+    TGeoVolume *itsV = gGeoManager->GetVolume(gm->GetITSVolPattern());
+    if (!itsV) {
+      printf("ITS volume %s is not in the geometry\n", gm->GetITSVolPattern());
+      return;
+    }
+    //
+    // Loop on all ITSV nodes, count Layer volumes by checking names
+    Int_t nNodes = itsV->GetNodes()->GetEntries();
+    Int_t numberOfLayers = 0;
+    for (Int_t j=0; j<nNodes; j++) {
+      if ( strstr(itsV->GetNodes()->At(j)->GetName(),
+                 gm->GetITSWrapVolPattern()) ) {
+
+       TGeoNode *wrapN = (TGeoNode*)(itsV->GetNodes()->At(j));
+       TGeoVolume *wrapV = wrapN->GetVolume();
+       Int_t nNodesWrap = wrapV->GetNodes()->GetEntries();
+
+       for (Int_t k=0; k<nNodesWrap; k++) {
+         if ( strstr(wrapV->GetNodes()->At(k)->GetName(),
+                     gm->GetITSLayerPattern()) ) {
+
+           Int_t l = numberOfLayers;
+           numberOfLayers++;
+       
+           char laynam[30];
+           snprintf(laynam, 30, "%s%d", gm->GetITSLayerPattern(), l);
+           TGeoVolume* volLy = gGeoManager->GetVolume(laynam);
+           printf("%s\t",volLy->GetName());
+
+           // Layers are Assemblies not Tubes, so:
+           // for rmax we assume the maximum extension of the assembly
+           // for rmin: for OB (no Turbo) we take the stave height and
+           // we subtract this from rmax; for IB (Turbo) we compute the
+           // position of the most internal corner
+           TGeoBBox *b = (TGeoBBox*)(volLy->GetShape());
+           TGeoBBox *s = (TGeoBBox*)(gGeoManager->GetVolume(Form("%s%d",gm->GetITSStavePattern(),l))->GetShape());
+           Double_t minr;
+
+           Double_t loc[3], mas[3];
+           if (l < 3) {
+             Double_t loc[3], mas[3];
+             loc[0] = s->GetDX();
+             loc[1] = s->GetDY();
+             loc[2] = 0;
+             volLy->GetNode(Form("%s%d_0",gm->GetITSStavePattern(),l))->GetMatrix()->LocalToMaster(loc,mas);
+             minr = TMath::Sqrt(mas[0]*mas[0] + mas[1]*mas[1]);
+             zPos[l+1] = 0;
+           } else {
+             minr = b->GetDX() - 2*s->GetDX();
+             TGeoBBox *c = (TGeoBBox*)(gGeoManager->GetVolume(Form("%s%d",gm->GetITSChipPattern(),l))->GetShape());
+             zPos[l+1] = c->GetDZ();
+           }
+           rmin[l+1] = minr;
+           rmax[l+1] = b->GetDX();
+           nlad[l+1] = gm->GetNStaves(l);
+
+           printf(" r = (%5.2lf , %5.2lf) ", rmin[l+1], rmax[l+1]);
+           printf(" z = +/- %5.2lf ", b->GetDZ());
+           printf(" #lad = %d \n", nlad[l+1]);
+
+         }
+       }
+      }
+    }
+  }
+
+}
+
+
+//______________________________________________________________________
+void DrawMaterialBudget_inPhi(const char *geofile = "geometry.root",
+                             Int_t lay = -1, Bool_t fromGDML = 0)
+{
+  //
+  // Function to compute material budget as a function of phi
+  // for layer lay (all layers and pipe if lay=-1). geofile is the
+  // geometry file name. If fromGDML is true, use a GDML file, else
+  // a standard geometry.root file.
+  // Cleaned up and simplified - M.Sitta 18 Nov 2014
+  //
+
+  Double_t rmin[8], rmax[8], zPos[8];
+  Int_t nlad[8];
+
+  ComputeGeometricalParameters(geofile, rmin, rmax, zPos, nlad, fromGDML);
+
+  Int_t n = 100000; //testtracks
+  // Test Get material ?
+  Double_t mparam[7]={0.,0.,0.,0.,0.};
+  Double_t point1[3],point2[3];
+
+  TH2F *xOvsPhi[8];
+  Int_t ls = 0, le = 8;
+  if (lay != -1) {
+    ls = lay+1;
+    le = lay+2;
+  }
+
+  for(Int_t i=ls; i<le; i++) {
+    Double_t phiMax;
+    if (i == 0) // PIPE doesn't have staves
+      phiMax = TMath::TwoPi()/4.;
+    else
+      phiMax = TMath::TwoPi()/nlad[i]*5.; // show approx 5 ladders
+    
+//    xOvsPhi[i] = new TH2F("", "", 100, 0, phiMax, 100, 0.2, 0.8);
+    if (i < 4)
+      xOvsPhi[i] = new TH2F("", "", 100, 0, phiMax, 100, 0.0, 1.2);
+    else
+      xOvsPhi[i] = new TH2F("", "", 100, 0, phiMax, 100, 0.2, 2.2);
+
+    for (Int_t it=0; it<n; it++) {
+      Double_t phi = phiMax*gRandom->Rndm();
+      Double_t z = zPos[i];
+      point1[0] = rmin[i]*TMath::Cos(phi);
+      point1[1] = rmin[i]*TMath::Sin(phi);
+      point1[2] = z;
+      point2[0] = rmax[i]*TMath::Cos(phi);
+      point2[1] = rmax[i]*TMath::Sin(phi);
+      point2[2] = z;
+      AliTracker::MeanMaterialBudget(point1,point2,mparam);
+      if (mparam[1] < 100)  // don't save fakes due to errors
+       xOvsPhi[i]->Fill(phi,mparam[1]*100); //fxOverX0Layer
+      if (!(it%10000)) cout << "layer" << i-1 << " : " << mparam[1]*100
+                           << " r=(" <<rmin[i] << "," << rmax[i] << ")"
+                           << " phi=" << phi <<endl;
+    }
+
+  }
+
+  if (lay==-1) {
+    TCanvas *c1 = new TCanvas("c1","Material Budget",700,800);
+    c1->Divide(2,4);
+    for(Int_t i=0; i<8; i++) {
+      c1->cd(i+1);
+      if (i>0) 
+       xOvsPhi[i]->SetTitle(Form("Layer %d",i-1));
+      else
+       xOvsPhi[i]->SetTitle(Form("Beampipe"));
+      xOvsPhi[i]->GetXaxis()->SetTitle("phi (rad)");
+      xOvsPhi[i]->GetYaxis()->SetTitle("X/X_{0} (%)");
+      xOvsPhi[i]->SetStats(0);
+      xOvsPhi[i]->Draw("col");
+    }
+  } else {
+    TCanvas *c1 = new TCanvas("c1","Material Budget");
+    Int_t i = lay+1;
+    xOvsPhi[i]->SetTitle(Form("Layer %d",lay));
+    xOvsPhi[i]->GetXaxis()->SetTitle("phi (rad)");
+    xOvsPhi[i]->GetYaxis()->SetTitle("X/X_{0} (%)");
+    xOvsPhi[i]->SetStats(1111111);
+    xOvsPhi[i]->Draw("col");
+  }
+
+}
+
+
+//______________________________________________________________________
+void DrawMaterialBudget_EtaVsPhi(const char *geofile = "geometry.root",
+                                Int_t lay = -1, Bool_t fromGDML = 0)
+{
+  //
+  // Function to compute material budget as a function of eta and phi
+  // in eta +/- 1 for layer lay (all layers and pipe if lay=-1).
+  // geofile is the geometry file name. If fromGDML is true, use a
+  // GDML file, else a standard geometry.root file.
+  // Cleaned up and simplified - M.Sitta 20 Nov 2014
+  //
+
+  Double_t rmin[10], rmax[10], zPos[10]; // zPos not used but needed for call
+  Int_t nlad[10];
+
+  ComputeGeometricalParameters(geofile, rmin, rmax, zPos, nlad, fromGDML);
+
+  rmin[8] = rmin[1];
+  rmax[8] = rmax[3];
+  rmin[9] = rmin[4];
+  rmax[9] = rmax[7];
+
+  Int_t n = 100000; //testtracks
+  // Test Get material ?
+  Double_t mparam[7]={0.,0.,0.,0.,0.};
+  Double_t point1[3],point2[3];
+
+  TH2F *xOverX0[10];
+  TH2F *nParticles[10];
+  Int_t ls = 0, le = 10;
+  if (lay != -1) {
+    ls = lay+1;
+    le = lay+2;
+  }
+
+  for(Int_t i=ls; i<le; i++) {
+    Double_t phiMax;
+    if (i == 0) // PIPE doesn't have staves
+      phiMax = TMath::TwoPi()/4.;
+    else if (i == 8) // Mean value of layer 0 to 2
+      phiMax = TMath::TwoPi()/(0.5*(nlad[1]+nlad[3]))*5.;
+    else if (i == 9) // Mean value of layer 3 to 6
+      phiMax = TMath::TwoPi()/(0.5*(nlad[4]+nlad[7]))*5.;
+    else
+      phiMax = TMath::TwoPi()/nlad[i]*5.; // show approx 5 ladders
+
+    xOverX0[i]    = new TH2F("", "", 100, -1., 1., 100, 0., phiMax);
+    nParticles[i] = new TH2F("", "", 100, -1., 1., 100, 0., phiMax);
+
+    for (Int_t it=0; it<n; it++) {
+      Double_t phi = phiMax*gRandom->Rndm();
+      Double_t eta = gRandom->Rndm()*2 - 1.; // +/- 1 eta
+      Double_t theta = TMath::ATan(TMath::Exp(-eta))*2;
+      point1[0] = rmin[i]*TMath::Cos(phi);
+      point1[1] = rmin[i]*TMath::Sin(phi);
+      point1[2] = rmin[i]/TMath::Tan(theta);
+      point2[0] = rmax[i]*TMath::Cos(phi);
+      point2[1] = rmax[i]*TMath::Sin(phi);
+      point2[2] = rmax[i]/TMath::Tan(theta);
+
+      AliTracker::MeanMaterialBudget(point1,point2,mparam);
+      if (mparam[1] < 100) { // don't save fakes due to errors
+       xOverX0[i]->Fill(eta,phi,mparam[1]*100); //fxOverX0Layer
+       nParticles[i]->Fill(eta,phi,1.);
+      }
+      if (!(it%10000)) cout << "layer" << i-1 << " : " << mparam[1]*100
+                           << " r=(" <<rmin[i] << "," << rmax[i] << ")"
+                           << " phi=" << phi << " theta=" << theta <<endl;
+    }
+
+    // normalization to number of particles
+    for (Int_t ix=1; ix<=nParticles[i]->GetNbinsX(); ix++) {
+      for (Int_t iy=1; iy<=nParticles[i]->GetNbinsY(); iy++) {
+       if (nParticles[i]->GetBinContent(ix,iy) > 0) 
+         xOverX0[i]->SetBinContent(ix,iy,xOverX0[i]->GetBinContent(ix,iy)/nParticles[i]->GetBinContent(ix,iy));
+      }
+    }
+  }
+  
+  if (lay==-1) {
+    TCanvas *c1 = new TCanvas("c1","Material Budget",750,900);
+    c1->Divide(2,5);
+    for(Int_t i=0; i<10; i++) {
+      c1->cd(i+1);
+      if (i==0) 
+       xOverX0[i]->SetTitle(Form("Beampipe"));
+      else if (i==8)
+       xOverX0[i]->SetTitle(Form("Layer 0 to 2"));
+      else if (i==9)
+       xOverX0[i]->SetTitle(Form("Layer 3 to 6"));
+      else
+       xOverX0[i]->SetTitle(Form("Layer %d",i-1));
+      
+      
+      xOverX0[i]->SetStats(0);
+      xOverX0[i]->GetXaxis()->SetTitle("pseudorapidity #eta ");
+      xOverX0[i]->GetYaxis()->SetTitle("#phi (rad)");
+      xOverX0[i]->GetZaxis()->SetTitle("X/X_{0} (%)");
+      xOverX0[i]->Draw("colz");
+    }
+  } else {
+    TCanvas *c1 = new TCanvas("c1","Material Budget");
+    Int_t i = lay+1;
+    xOverX0[i]->SetTitle(Form("Layer %d",lay));
+    xOverX0[i]->SetStats(0);
+    xOverX0[i]->GetXaxis()->SetTitle("pseudorapidity #eta ");
+    xOverX0[i]->GetYaxis()->SetTitle("#phi (rad)");
+    xOverX0[i]->GetZaxis()->SetTitle("X/X_{0} (%)");
+    xOverX0[i]->Draw("colz");
+  }
+
+}
diff --git a/ITS/UPGRADE/testITSUv1/MakeMatBudPlots.C b/ITS/UPGRADE/testITSUv1/MakeMatBudPlots.C
new file mode 100644 (file)
index 0000000..5471d22
--- /dev/null
@@ -0,0 +1,46 @@
+void MakeMatBudPlots(int layer=0)
+{
+  // Simple interface to GetMaterialBudget.C macro to plot the
+  // material budget (for each material separately) for a given layer
+  // M.S. 01 Jul 2014
+
+  const int kNLr = 7;
+  const int kNLrInner = 3;
+
+  gSystem->Load("libITSUpgradeBase");
+  gSystem->Load("libITSUpgradeSim");
+  gROOT->LoadMacro("GetMaterialBudget.C");
+
+  enum {kRmn,kRmd,kRmx,kNModPerStave,kPhi0,kNStave,kNPar};
+  // Radii are from last TDR (ALICE-TDR-017.pdf Tab. 1.1, rMid is mean value)
+  // TO BE KEPT IN SYNC WITH CreateITSUv1.C MACRO!!!
+  const double tdr5dat[kNLr][kNPar] = { 
+    {2.24, 2.34, 2.67,  9., 16.37, 12}, // for each inner layer: rMin,rMid,rMax,NChip/Stave, phi0, nStaves
+    {3.01, 3.15, 3.46,  9., 12.03, 16},
+    {3.78, 3.93, 4.21,  9., 10.02, 20},
+    {-1,  19.6 ,   -1,  4.,  0.  , 24},  // for others: -, rMid, -, NMod/HStave, phi0, nStaves
+    {-1,  24.55, -1,    4.,  0.  , 30},
+    {-1,  34.39, -1,    7.,  0.  , 42},
+    {-1,  39.34, -1,    7.,  0.  , 48} 
+  };
+
+  if (layer < 0 || layer >= kNLr) {
+    printf("Wrong layer number %d - giving up\n",layer);
+    return;
+  }
+
+  double rmin, rmax;
+
+  if (layer < kNLrInner) { // Inner layers
+    rmin = tdr5dat[layer][kRmn] - 0.1;
+    rmax = tdr5dat[layer][kRmx] + 0.1;
+  } else { // Outer layers
+    rmin = tdr5dat[layer][kRmd] - 0.175;
+    rmax = tdr5dat[layer][kRmd] + 0.35;
+  }
+
+  printf("Drawing material budget for layer %d (Rmin = %f Rmax = %f\n",layer,
+        rmin,rmax);
+  DrawMaterialBudget_Splitted(layer, rmin, rmax) ;
+
+}
diff --git a/ITS/UPGRADE/testITSUv1/MaterialBudget_README b/ITS/UPGRADE/testITSUv1/MaterialBudget_README
new file mode 100644 (file)
index 0000000..460d647
--- /dev/null
@@ -0,0 +1,120 @@
+Simple instructions on how to produce Material Budget plots
+===========================================================
+
+Here there are some instruction on how to produce Material Budget
+plots showing separately the contribution of each component or
+material (such as, e.g., Figs. 4.3 in the TDR):
+
+
+1) all needed files are in $ALICE_ROOT/ITS/UPGRADE ; in particular
+use the Config.C and CreateITSUv1.C in there
+
+
+2) geometry files with each component are needed; to produce them
+it is enough to run a simulation changing the value of
+
+        const int kBuildLevel
+
+in CreateITSUv1.C ; kBuildLevel = 0 corresponds to the complete
+description with all materials (and should be used for any real
+simulation). Geometry files obtained with values of kBuildLevel from
+0 to 5 (for the Inner Barrel) or 6 (for the Outer Barrel) are needed
+(for the IB the build levels 5 and 6 are identical).
+
+
+3) to easily produce all needed file in one go, the script runMatBud.sh
+can be used. Simply execute it to produce all geometry.root files
+needed for further analysis. The file CreateITSUv1.C_template is
+required in the same directory where the script is executed: it sets the
+proper value of kBuildLevel and saves it as CreateITSUv1.C to be used in
+the simulation.
+
+The runMatBud.sh script accepts an optional parameter, the maximal
+kBuildLevel value to be used; the default is 6. So
+
+        ./runMatBud.sh
+
+produces seven geometry files (named from geometry_0.root to
+geometry_6.root), while for instance
+
+        ./runMatBud.sh 5
+
+produces six geometry files (named from geometry_0.root to geometry_5.root)
+
+
+4) Once the geometry files are produced, the macro GetMaterialBudget.C
+can be used to produce the plots. After starting an aliroot session,
+first of all the needed libraries and the macro are to be loaded
+
+        aliroot
+        root [0] gSystem->Load("libITSUpgradeBase")
+        root [1] gSystem->Load("libITSUpgradeSim")
+        root [2] .L GetMaterialBudget.C
+
+then the function DrawMaterialBudget_Splitted can be called to
+produce a plot similar to those of the TDR, i.e. the X/X0 fraction
+as a function of the azimuthal angle phi for a given layer. The
+function takes three arguments: the layer number, the minimum radius
+and the maximum radius, for instance
+
+        root [3] DrawMaterialBudget_Splitted(3,19.4,19.9)
+
+the macro takes care automatically of the different number of build
+levels for IB and OB. The plot is also saved in a file named
+Material-details.pdf
+
+Other useful functions in the same macro are
+i) PrintMaterialDefs to print a list of all ITS materials and their
+definition; it takes the geometry.root file as optional argument
+
+        root [4] PrintMaterialDefs()
+
+ii) DrawMaterialBudget_FromTo to plot a 2D histogram with the material
+budget for a given layer vs eta and phi. Warning: the function assumes
+that the geometry file is called "geometry.root"!
+
+The function takes four arguments: the first three are the layer number
+and the minimum and maximum radius; if the optional fourth argument
+is false (kFALSE or 0) also a 1D histogram is produced with the budget
+as a function of phi (the default is true, so only the 2D plot is shown)
+
+        root [5] DrawMaterialBudget_FromTo(4,24.3,24.8,0)
+
+The plot is also saved in a file named either Material-2D.pdf or
+Material-1D.pdf
+
+iii) DrawMaterialBudget_inPhi to produce 2D plots of material budget
+vs phi and fraction of X/X0. The macro takes three optional arguments:
+the geometry file (defaults to "geometry.root"), the layer number and
+a flag to use a Root file (default) or a GDML file (if true)
+
+        root [6] DrawMaterialBudget_inPhi()
+
+if no layer is specified, all layers and the beam pipe are shown
+
+iv) DrawMaterialBudget_EtaVsPhi to produce a 2D plot of material
+budget as a function of eta (in the range +/- 1) and phi. The macro
+takes three optional arguments: the geometry file (defaults to
+"geometry.root"), the layer number and a flag to use a Root file
+(default) or a GDML file (if true)
+
+        root [7] DrawMaterialBudget_inPhi()
+
+if no layer is specified, all layers and the beam pipe are shown
+along with the total material for the Inner and Outer Barrels
+separately.
+
+
+5) To avoid the need of loading by hand the libraries and the macro
+and the necessity of entering the correct radii, the same plot of
+the TDR produced by the function DrawMaterialBudget_Splitted can be
+obtained with the macro MakeMatBudPlots : it takes only one argument,
+the layer number, and takes care of the housekeeping needed to get
+the plot
+
+        root [7] .x MakeMatBudPlots.C(2)
+
+
+-----------------------------------------------------------
+For questions, comments or suggestions: Mario.Sitta@cern.ch
+
diff --git a/ITS/UPGRADE/testITSUv1/runMatBud.sh b/ITS/UPGRADE/testITSUv1/runMatBud.sh
new file mode 100755 (executable)
index 0000000..02428fd
--- /dev/null
@@ -0,0 +1,24 @@
+#!/bin/bash
+
+# This simple procedure produces all geometry_*root files
+# needed for an analysis of the material budget using the
+# macro GetMaterialBudget.C and MakeMatBudPlots.C
+# It takes one optional parameter, the maximum  number of build levels.
+# It only requires CreateITSUv1.C_template as input
+# M.Sitta 11 Nov 2014
+
+if [ "x$1" != "x" ]; then
+  MAXBUILDLEVEL=$1
+else
+  MAXBUILDLEVEL=6
+fi
+
+rm -Rf geometry_[0-5].root
+
+for i in `seq 0 $MAXBUILDLEVEL`; do
+  sed "s/BUILDLEVEL/$i/" CreateITSUv1.C_template >CreateITSUv1.C
+  aliroot -b -q sim.C\(1\)
+  mv geometry.root geometry_$i.root
+done
+
+ls -l geometry_*root