Example 5 Implementation of the FLUKA simulation of the Pb-Scintillating fiber e.m. calorimeter of the KLOE Experiment -------------------------------------------------------- UPDATED===26.10.2010.14.10.46 TITLE===Implementation of the FLUKA simulation of the Pb-Scintillating fiber e.m. calorimeter of the KLOE Experiment TYPE===database -------------------------------------------------------- --------------------------------------------------------

Giuseppe Battistoni; INFN, Milano
Jaroslaw Zdebik; UJ, Cracow

--------------------------------------------------------

J. Zdebik

Abstract:

In this example, I describe the process and results of implementation geometry of the KLOE Barrel Calorimeter with realistic material composition into Fluka package. All my work is based on implemented one calorimeter module with rectangular shape (work was done by: "G. Battistoni, B. Di Micco, A. Ferrari, A. Passeri, V. Patera" ) inside Fluka framework. I describe shortly the development of this module in order to achieve the trapezoid shape (which exists in the realistic world) and finally, the implementation of the complete structure of twenty four modules geometry of the Kloe Barrel Calorimeter (work was done by: J. Zdebik with "B. Di Micco suggestions" )



Introduction. ::::::::::::

  In this example I will present the the most important steps and code fragments which were necessary to perform in order to achieve the complete simulation of the KLOE Barrel Calorimeter geometry with FLUKA, also implementation of the realistic calorimeter material composition will be discussed. I will focus only on the code related with the FLUKA software.


The KLOE calorimeter. ::::::::::::

  The KLOE calorimeter consists of three parts: one barrel and two endcaps. Each part is built out of modules with a spaghetti structure made of lead layers with about 1,2 mm thickness (200 layers per 1 module) which are filled with scintillating fibers of 1 mm diameter. The barrel (Fig. 1) has a cilindrical shape and is divided into 24 modules with trapezoid shape. Each module is 23 cm thick.

KLOE barrel calorimeter
Figure 1: A photo of the KLOE calorimeter [1].
The 24 modules of the barrel and the inner plane of the endcap are shown.

The volume of the calorimeter consists of 50% fiber, 40% lead and 10% glue. The section of one of the 24 modules and the details of the lead+fiber structure can be aprreciated in Fig. 2. The detailed structure of the module has been simulated as described in the following section.

trapes_structure fibers_structure
Figure 2: Single module section (left) and details of the fibers structure in the calorimeter detector (right).

Material description. ::::::::::::

  The calorimeter is built out of lead, plexi and glue [2] which are composed as follows:

The hardener has been simulated as a mixture of several materials which composition is given in Table 1.

Table 1: Composition of hardener mixture [2].

Compound formula fraction
Polyoxypropylediamine C7H20O3 90%
Triethanolamine C6H15O3 7%
Aminoethylpiperazine C6H15N3 1.5%
Diethylenediamine C4H10N2 1.5%

Single module geometry implementation. ::::::::::::

  The module geometry was designed trying to profit of the big degree of symmetry of the calorimeter structure. The first version of this geometry was realized in rectangular shape [Implementation to FLUKA package was done by: G. Battistoni, B. Di Micco, A. Ferrari, A. Passeri and V. Patera ] with the material description given in the previous section.
FLUKA allows to define a symmetric geometric structure as recursive isometric transformations of a base structure. The base structure (referred to as base cell in the following) has to be defined in its shape and material composition, while the transformed regions (referred to as replicas in the following) are defined as void boxes whose boundary surfaces can be overlapped to the boundaries of the base cell. At the tracking level each time a particle enters in a replica region a transformation is performed in the base cell. The particle is propagated in the base cell and transformed back to the replica when it reachs the boundary of the base cell. These transformations, forward and backward, have to be provided by the user in the lattic.f routine whose code is reported in source code area .
In the reference system shown in Fig. 3 (the blue axis denote the frame which was used for simulations with FLUKA, the red coordinate system presents the coordinate frame of the KLOE detector) the base cell has dimensions of 52 cm x 1.2 mm x 430 cm, it consists of a lead block filled with 385 scintillating fibers and glue cylinders (Fig. 4).

base_cell
Figure 3: The coordinate systems used in simulations.
base_cell
Figure 4: Details of the implementation of the calorimeter layers [3].

This cell is replicated 199 times to build a total number of 200 layers forming the calorimeter module (see left panel in Fig. 5).
trapezoid_shape rectangular_shape
Figure 5: Cross section of the EmC module with rectangular and trapezoid geometry.
Visualization by the FLAIR program [4] using as an input a geometry setup files of FLUKA.

Each plane is alternatively shifted of 0.675 mm in order to reproduce the real fiber configuration (Fig. 6).
base_cell
Figure 6: The KLOE fibers schematic view.

  In order to implement the whole barrel calorimeter into FLUKA we extended the rectangular module to the realistic trapezoid shape. This has been achieved by declaration of a new base cell consisting of a scintillating fiber and glue cylinder both inserted into a small block of lead. The length of the cell is equal to the length of the module (4.3 m) and its cross section corresponds to the square of 1.2 mm, where the fiber diameter is equal to 1 mm and width of the glue cylinder equals to 0.1 mm.
Then this new base cell was replicated about 4000 times to build two triangular sections on the left and right side of the module (see Fig. 7). A lattice transformation from the 4000 replica to the base cell provides the needed translation.
base_cell
Figure 7: Structure of the single module with trapezoid shape.

  A complete visualisation of the new trapezoid geometry [Visualization is based on the FLUKA sources and prepared using a FLAIR (FLUKA advanced visualization geometry tool [4])] is shown in Fig. 5, whereas the scatter-plot of single energy deposits from homogenously distribution photons hittting the module are shown in Fig. 8: a trapezoid structure which was built with realistic energy deposits from the scintillating fibers can be easily recognized.
energy deposits
Figure 8: Distribution of energy deposits in scintillating fibers
of a single unit of the electromagnetic calorimeter.

  As a next step, it would be natural to replicate 23 times the defined geometry for one module and to build a full barrel. This task was unfortunately impossible to realise with the using by us version of FLUKA [FLUKA 2008.3 for GNU linux operating system] because this version doesn't permit to replicate the region which had been replicated before [5]. It is due to the fact that the lattice replication on the second and higher levels is not implemented yet. This is the reason why geometry was built in another, unfortunatelly much more complicated way.

The barrel geometry implementation. ::::::::::::

  The barrel was first defined as 24 empty volumes (containers). In one of these volumes we implemented the base cell which describes one layer of lead, fibers and glue. This cell was replicated 199 times in this volume, to build a rectangular section part of the main trapezoid module (Fig. 9, left), the same idea was used to build a single module with rectangular shape, and this base cell was replicated 200 times in each of the remaining 23 modules (see Fig. 10). A program code with lattice transformation that performs a rotation and 200 translations for each module can be examined in appendix A. The remaining two triangle areas in the first module are filled with the real structure of lead, fibers and glue cylinder. These two regions were then replicated at the corresponding positions in the remaining modules, through a single rotation for each module (Fig. 9, left)
  The left panel of Fig. 9 presents details of the edge area between the base module and the first replicated area on its right side. In the right panel of Fig. 9 the energy deposits in the scintillating fibers at the edge of modules are presented. The figure was obtained for uniformly distributed photons.

energy deposits energy deposits
Figure 9: Particular of the base cell and details of the implementation on the edge of two modules.

Fig. 10 shows a visualisation of the whole barrel calorimeter geometry. Using this geometry we are able to study in detail energy deposits also at the edges of the modules.
energy deposits
Figure 10: Visualization of the 24 trapezoid modules of the barrel calorimeter with FLAIR.

Finally in Fig. 11, as an example of the proper implementation of the geometry we show energy deposits in scintillating fibers in the whole barrel calorimeter.
geometry details geometry details
Figure 11: Energy depositions in fibers in the KLOE barrel calorimeter.


Light output simulation. ::::::::::::

  The energy deposits are computed by FLUKA [6] taking into account the Birks effect [7] (see equation 1), that is the saturation of the light output (L) of a scintillating material when the energy release is high [6]. For high densities of energy deposition, due to the quenching interactions between the excited molecules along the path of incident particles, the light output is not changing linearly with the energy deposition but instead it can be described as [8]:

Birk [1]
This law describes the light output of (organic) scintillators [8]. The luminescence at low specific ionisation density depends on the particle type and the material, and the parameters c1 and c2 are equal to [3]:
c1 = 0.013 cm2MeV-1
c2 = 9.6 x 10-6 cm3MeV-2
A proper attenuation length parameter for fibers material in order to calculate a light intensity at each side of the module was taken into account. The attenuation function of fibers was described with the following formula:
attenuation [2]
where y is distance between the place of deposited energy and the photocathode. The attenuation factor B is the ratio between detected and generated light signal which changes as a function of the distance (on y axis) from the generation point to the photocatode. The module consists of two types of scintillating fibers: "Kuraray" and "Pol.hi.tech" [2]. The first are implemented until a depth of 12 cm and the remaining part of the module is built with latter type of fibers. The value of the parameters for "Kuraray" and "Pol.hi.tech" are shown in Table 2.

Table 2: Fibers parameters.

A lambda_1 [cm] lambda_2 [cm]
Kuraray 0.35 50 430
Pol.hi.tech 0.35 50 330

Some parts of source code. ::::::::::::

A) Building one module with trapezoid shape.

1. WORK WITH 1 module with RECTANGULAR shape.

1.1 Almost all parts of this geometry are the same (as were develped
  for single module with rectangular shape), structure of lattices, and base cell in
  rectangular part, (more informations you can find in [2] and [9]
  positions in bibliography). 

Modification was based on shifting a lattices 1,3,5,7,.....(odd numbers) on x axis on
0.135/2.0 cm in geometry file in order to produce the following
structure of fibers: 

o o o o o
 o o o o o
o o o o o

z
/\
|
|
-------> x

2. Description:  

2.1 Implementation a new base cell with:

lead block
  RPP  774 0.520350000000000E+02 0.521700000000000E+02-0.215000000000000E+03
           0.215000000000000E+03 0.105376086559900E+01 0.117067429511000E+01

one glue cylinder                                                              
* bdjj glue 
  YCC  775 1.112217580350000E+00 0.521025000000000E+02 0.550000000000000E-01                                                                          
one fiber cylinder
*  bdjj  fiber                                                             
  YCC  776 1.112217580350000E+00 0.521025000000000E+02 0.500000000000000E-01

These cylinders were located in the triangle areas:
22.....3 2 1                1 2 3 4 ..... 22
*   -------------------------------                                               
*    ooooo|                 |ooooo                                         
*      oooo|                 |oooo                                           
*      ooo|                 |ooo                                           
*        oo|                 |oo                                            
*        o|                 |X                                            
*          |                 |
*          -------------------                                               
*                                                                          
*                                                                          
* X - base cell which was replicated          
* o - new lattices = 199 - 9 = 190                                         
*                                                                          

New lattices are situated in 22 columns on the left side of old rectangular part
and the same 22 columns for the right side (one box = a new base cell).

****************************************************************************
*             new lattices on the right and left side			   *
****************************************************************************
*
*     777  lattice     10      lattice right side column 1
  RPP  777 0.519675000000000E+02 0.521025000000000E+02-0.215000000000000E+03
           0.215000000000000E+03 0.117067429511000E+01 0.128758772462100E+01
.
.
.
.
.
.
*    5020  lattice     198     lattice left side column 22
  RPP 5020-0.297750000000000E+01-0.284250000000000E+01-0.215000000000000E+03
           0.215000000000000E+03 0.231503990431781E+02 0.232673124726891E+02
*    5021  lattice     199     lattice left side column 22
  RPP 5021-0.291000000000000E+01-0.277500000000000E+01-0.215000000000000E+03
           0.215000000000000E+03 0.232673124726891E+02 0.233842259022001E+02
*
*
* Each lattice is declared as a region.
*

***************************************************************************
****** Old 199 lattices in rectangular shape - modificated - new coordinates.
*
*
*     5022  lattice     1
  RPP 5022 0.600000000000000E-01 0.520350000000000E+02-0.215000000000000E+03
           0.215000000000000E+03 0.118453429511000E+00 0.235366859022000E+00
.
.
.
.
.
*     5220  lattice     199
  RPP 5220 0.600000000000000E-01 0.520350000000000E+02-0.215000000000000E+03
           0.215000000000000E+03 0.232673124726891E+02 0.233842259022001E+02
*
*

2.2. Input file:

* Declared all new numbers of the regions with material
* We are using GLOBAL card:
   GLOBAL    5000.0

*------------ new modyficated code by JJ ---------------------------*
*
*----- old fibers regions
*
ASSIGNMAT   23.0    6.0       390.0               1.0
*
*----- old glue region
*
ASSIGNMAT   31.0    391.0                         1.0
*
*----- old vacuum region
*
ASSIGNMAT   2.0     392.0
*
*----- old lattices
*
ASSIGNMAT   24.     393.0     591.0               1.0
*
*--- new parts for trapesoidal areas ----------------------------------*
*
ASSIGNMAT   24.     592.0     4836.0              1.0
* lead block
ASSIGNMAT   22.0    3.0                           1.0
*  775 - glue
ASSIGNMAT   31.0    5.0                           1.0
*  776 - fiber
ASSIGNMAT   23.0    4.0                           1.0
*
*--------------------------------------------------------------*

For a new lattice areas was not needing to declare a materials.

2.3. mgdraw.f routine:

- Declarations a lfibra variable tru for all regions with fibers,
- One fiber is contained in 1 region.

2.4 lattic.f routine:

- For the old part is the same,
- For the new part is based on using IF conditions. 

------------- For example:

*  First column on the right (899  + 189  = 1088)
*  not + 190 (because the first is a base cell)
*
      IF (IRLTGG.GT.(898).AND.IRLTGG.LT.(1089)) THEN
        IF(MOD((IRLTGG),2).EQ.0) THEN
            SB (1) = XB (1)
        ELSE
            SB (1) = XB (1) + 0.135D+00/TWOTWO
        ENDIF
          SB (2) = XB (2)
          SB (3) = XB (3)-1.169134295110000D-01*(1+IRLTGG-899)
          UB (1) = WB (1)
          UB (2) = WB (2)
          UB (3) = WB (3)

      END IF

*
*  Second column on the right (1089 + 181  = 1270)
*
      IF (IRLTGG.GT.(1088).AND.IRLTGG.LT.(1271)) THEN
       IF(MOD((IRLTGG),2).EQ.0) THEN
        SB (1) = XB (1) - 0.135D+00
       ELSE
        SB (1) = XB (1) + 0.135D+00/TWOTWO - 0.135D+00
       ENDIF
        SB (2) = XB (2)
        SB (3) = XB (3)-1.169134295110000D-01 * (9+IRLTGG - 1089)
        UB (1) = WB (1)
        UB (2) = WB (2)
        UB (3) = WB (3)

      END IF

------------------- 

2.5 Some details:

TOTAL NUMBER of regions: 4836 
TOTAL NUMBER of  bodies: 5220

NUMBER of lines in geometry file: 21230
B) Lattice replication for 24 modules.
 
Below is presented the source code of the lattice transformation
routine, which build in fly the whole barrel calorimeter geometry, 
this is the heart of the whole simulation which allows to build barrel calorimeter geometry in FLUKA simulation.

*=== lattic ===========================================================*
*
      SUBROUTINE LATTIC (XB, WB, DIST, SB, UB, IR, IRLTGG, IFLAG)

      INCLUDE '(DBLPRC)'
      INCLUDE '(DIMPAR)'
      INCLUDE '(IOUNIT)'
      INCLUDE '(RTGMMV)'

c      INCLUDE '(GLTLOC)'

*----------------------------------------------------------------------*
*     LATTIC: user written routine which must return the tracking point
*     and direction ( SB, UB ) corresponding to region number IR, cell 
*     number IRLTGG and real position/direction XB, WB                 
*                                                                      
*     Created on 16 December 1993  by    Alfredo Ferrari & Paola Sala 
*                                                   Infn - Milan       
*                                                                      
*     Last change on 26-may-98     by    Alfredo Ferrari               
*----------------------------------------------------------------------*

      DIMENSION XB (3), WB (3), SB (3), UB (3), UN (3)

      DOUBLE PRECISION UNN(3),UBB(3)
      include 'kloecal.inc'

      LOGICAL LFIRST
      DATA LFIRST /.TRUE./

      IF (LFIRST) THEN
       LFIRST = .FALSE.
      ENDIF

c  rewritting to variable which we can modificate
c  X,W - where particle hit
c  S,U - response in the lattice cell

      KMOD=0

      SB(1)=XB(1)
      SB(2)=XB(2)
      SB(3)=XB(3)

      UB(1)=WB(1)
      UB(2)=WB(2)
      UB(3)=WB(3)

      IF(IRLTGG.lt.1.or.IRLTGG.gt.4822) THEN
       WRITE (LUNOUT,*)' NON-EXISTENT LATTICE: ',IRLTGG

c rotation of 23 trapesoids to tha base trapesoid

      ELSEIF (IRLTGG.GE.4800.and.IRLTGG.LE.4822) THEN
       XRT    =  SB(1)-XO
       YRT    =  SB(2)-YO
       ZRT    =  SB(3)-ZO
       THETA  =  -TWOTWO*PIPIPI/24.0D0*(IRLTGG-4799)

       CALL ROTATE(SB(1),SB(2),SB(3),XRT,YRT,ZRT,THETA)
        SB(1)=SB(1)+XO
        SB(2)=SB(2)+YO
        SB(3)=SB(3)+ZO
       CALL ROTATE(UBB(1),UBB(2),UBB(3),UB(1),UB(2),UB(3),THETA)
        DO i=1,3
         UB(i) = UBB(i)
        ENDDO

c transformation from 199 lattices into trapes 1

      ELSEIF (IRLTGG.GE.1.AND.IRLTGG.LE.199) THEN
        IF(MOD((IRLTGG),2).EQ.0) THEN
         SB (1) = SB (1)
        ELSE
         SB (1) = SB (1) - 0.135D+00/TWOTWO
        ENDIF

         SB (3) = SB (3) - 1.169134295110000D-01 * IRLTGG
c
c --------------------------------------------------
c transformation of the 200 layers for 2 - 24 trapesoids objects
c
* in which trapesoid I am
* KMOD from 2 to the 24, but rotation for 1 - 23
* KLAYER - layer number in each trapesoid 1-200
c

      ELSEIF (IRLTGG.GE.200.AND.IRLTGG.LE.4799) THEN
       KMOD   = (IRLTGG - 200)/200 + 1 + 1
       THETA  = -TWOTWO*PIPIPI/24.0D0*(KMOD - 1)
        KLAYER = IRLTGG - 200 * (KMOD-1) + 1
          XRT=SB(1)-XO
          YRT=SB(2)-YO
          ZRT=SB(3)-ZO
         CALL ROTATE(SB(1),SB(2),SB(3),XRT,YRT,ZRT,THETA)
          SB(1)=SB(1)+XO
          SB(2)=SB(2)+YO
          SB(3)=SB(3)+ZO
         CALL ROTATE(UBB(1),UBB(2),UBB(3),UB(1),UB(2),UB(3),THETA)
          DO i=1,3
           UB(i) = UBB(i)
          ENDDO

        IF(MOD((KLAYER),2).EQ.0) THEN
          SB (1) = SB (1) - 0.135D+00/TWOTWO
       ENDIF

*
* first layer without transformation because it is on the base cell level
*
*
        SB (3) = SB (3) - 1.169134295110000D-01 * (KLAYER-1)
     ENDIF


      RETURN

*=== lantor =======================================================*

      ENTRY LATNOR ( UN, IRLTNO )
*
*--------------------------------------------------------------------*
*   LATtice cell NORmal transformation:                              
*                                                                    
*   Input variables:  un(i) = normal components in the tracking re-  
*                               ference system                       
*                     irltno = present lattice cell \#               
*   Output variables:  un(i) = normal components in the problem re-  
*                                ference system                      
*--------------------------------------------------------------------*
*
c rotation of the 23 trapesoids

      IF (IRLTNO.GE.4800.and.IRLTNO.LE.4822) THEN
       THETA=+TWOTWO*PIPIPI/24.0D0*(IRLTNO-4799)

        CALL ROTATE(UNN(1),UNN(2),UNN(3),UN(1),UN(2),UN(3),THETA)
         DO I=1,3
          UN(I)=UNN(I)
         ENDDO

c for the 199 lattices in trapesoid-1

      ELSEIF (IRLTNO.GE.1.AND.IRLTNO.LE.199) THEN

c for the 1 - 200 layers in trapesoids number 2 - 24 

      ELSEIF (IRLTNO.GE.200.and.IRLTNO.LE.4799) THEN
       KMOD   = (IRLTNO - 200)/200 + 1 + 1
       THETA  = +TWOTWO*PIPIPI/24.0D0*(KMOD - 1)

        CALL ROTATE(UNN(1),UNN(2),UNN(3),UN(1),UN(2),UN(3),THETA)
         DO I=1,3
          UN(I)=UNN(I)
         ENDDO
      ENDIF

      RETURN
*=== End of subroutine lattic =======================================*

      END

*--------------------------------------------------------------------

      SUBROUTINE ROTATE(XOUT,YOUT,ZOUT,XIN,YIN,ZIN,THETA)
        DOUBLE PRECISION XIN,YIN,ZIN
        DOUBLE PRECISION XOUT,YOUT,ZOUT
        DOUBLE PRECISION THETA

        XOUT = COS(THETA)*XIN+SIN(THETA)*ZIN
        YOUT = YIN
        ZOUT = -SIN(THETA)*XIN+COS(THETA)*ZIN
      END


c===================================================================
c===================================================================

Conclusions. ::::::::::::
  We have built a full simulation of the Barrel Calorimeter of the KLOE detector with FLUKA Monte Carlo package. This implementation of the geometry was used for simulations of the physical energy response of the calorimeter. The Monte Carlo cooperates with vertex generator which reproduces kinematics of the physical reactions.
  The prepared program with geometry of the 24 modules of the barrel calorimeter enables fast simulations of a choosen reaction channel and facilitates studies of the detector response to different reactions separately.
  Results achieved with prepared simulation like: estimation influence of merging and splitting of clusters for reconstruction possibilities of the particles on the whole barrel calorimeter with a present clustering algorithm or investigation influence of inreasing quantum efficiency of photomultipliers on mentioned effects constituted a main results for my diploma thesis [9]. Bibliography. ::::::::::::
1. http://www.lnf.infn.it/kloe/
2. B. Di Micco et al., Nuclear Physics B (Proc. Suppl.) 172 (2007) 243-245
3. A. Ferrari et al., Status of the KLOE EmC Simulation with FLUKA , KLOE-2 meeting, LNF July 25th, (2006)
http://www.lnf.infn.it/kloe/kloe2/calorimeter/20060725_ferrari.ppt
4. FLAIR project homepage, http://www.fluka.org/flair/index.html
5. FLUKA organization homepage, http://www.fluka.org
6. M. Anelli et al., Nucl. Instr. & Meth. in Phys. Res. A 581 (2007) 368-372
7. J.B. Birks, The Theory and Practice of Scintillation Counting , Pergamon Press, London, (1964)
8. R. Leo, Techniques for Nuclear and Particle Physics Experiments , Springer-Verlag, Berlin Heidelberg, New York, London, Paris, Tokyo, (1987)
9. J. Zdebik, Merging and splitting of clusters in the electromagnetic calorimeter of the KLOE detector, LNF-08/19(Thesis) - July 18, (2008); JINST 2008_JINST_TH_004 (2008); arXiv:0811.1377v1 (2008), LAP Lambert (2010), in preparation