Skip to content
Snippets Groups Projects
ce_ecal_placement.py 7.62 KiB
Newer Older
  • Learn to ignore specific revisions
  • '''
        A script to calcualte placement of ecal endcap modules
        lxml is not included in container, get it by simply typing 'pip install lxml'
        Author: Chao Peng (ANL)
        Date: 06/17/2021
    '''
    
    import numpy as np
    import argparse
    from lxml import etree as ET
    
    
    # start pos, total number of blocks
    CRYSTAL_ALIGNMENT = [
        (8, 17),
        (8, 17),
        (8, 17),
        (7, 18),
        (7, 18),
        (6, 19),
        (4, 19),
        (1, 22),
        (1, 22),
        (1, 22),
        (1, 22),
        (1, 22),
        (1, 20),
        (1, 20),
        (1, 18),
        (1, 18),
        (1, 16),
        (1, 16),
        (1, 14),
        (1, 14),
        (1, 12),
        (1, 12),
        (1, 6),
        (1, 6),
    ]
    
    GLASS_ALIGNMENT = [
        (13, 11),
        (13, 11),
        (13, 11),
        (12, 12),
        (12, 12),
        (12, 12),
        (11, 12),
        (10, 13),
        (9, 13),
        (8, 14),
        (7, 14),
        (4, 16),
        (1, 19),
        (1, 18),
        (1, 18),
        (1, 17),
        (1, 17),
        (1, 15),
        (1, 13),
        (1, 11),
        (1, 10),
        (1, 8),
        (1, 6),
    ]
    
    # calculate positions of modules with a quad-alignment and module size
    def individual_placement(alignment, module_x=20.5, module_y=20.5):
        placements = []
        for row, (start, num) in enumerate(alignment):
            for col in np.arange(start - 1, start + num - 1):
                placements.append(((col + 0.5)*module_y, (row + 0.5)*module_x))
        placements = np.asarray(placements)
        return np.vstack((placements,
                np.vstack((placements.T[0]*-1., placements.T[1])).T,
                np.vstack((placements.T[0], placements.T[1]*-1.)).T,
                np.vstack((placements.T[0]*-1., placements.T[1]*-1.)).T))
    
    
    
    if __name__ == '__main__':
        parser = argparse.ArgumentParser()
        parser.add_argument('-s', '--save', default='compact/ce_ecal_crystal_glass.xml',
                help='path to save compact file.')
        args = parser.parse_args()
    
        data = ET.Element('lccdd')
        defines = ET.SubElement(data, 'define')
    
        # constants: name, value
        consts = [
            ('CrystalModule_sx', '20.0*mm'),
            ('CrystalModule_sy', '20.0*mm'),
            ('CrystalModule_sz', '20.0*cm'),
            ('CrystalModule_wrap', '0.5*mm'),
            ('CrystalModule_dz', '-8*cm'),
            ('GlassModule_sx', '40.0*mm'),
            ('GlassModule_sy', '40.0*mm'),
            ('GlassModule_sz', '40.0*cm'),
            ('GlassModule_dz', '0.0*cm'),
            ('GlassModule_wrap', '1.0*mm'),
            ('EcalEndcapN_z0', '-EcalEndcapN_zmin-max(CrystalModule_sz,GlassModule_sz)/2.'),
        ]
    
        for name, value in consts:
            constant = ET.SubElement(defines, 'constant')
            constant.set('name', name)
            constant.set('value', value)
    
        # this value will be used multiple times, so define it here
        readout_name = 'EcalEndcapNHits'
    
        # detector and its dimension/position/rotation
        dets = ET.SubElement(data, 'detectors')
        cmt = ET.SubElement(dets, 'comment')
        cmt.text = ' Backwards Endcap EM Calorimeter, placements generated by script '
    
        det = ET.SubElement(dets, 'detector')
        det.set('id', 'ECalEndcapN_ID')
        det.set('name', 'EcalEndcapN')
        det.set('type', 'HomogeneousCalorimeter')
        det.set('readout', readout_name)
    
        pos = ET.SubElement(det, 'position')
        pos.set('x', '0')
        pos.set('y', '0')
        pos.set('z', 'EcalEndcapN_z0')
    
        rot = ET.SubElement(det, 'rotation')
        rot.set('x', '0')
        rot.set('y', '0')
        rot.set('z', '0')
    
        # placements of modules
        plm = ET.SubElement(det, 'placements')
    
        # crystals
        # crystal = ET.SubElement(plm, 'individuals')
        crystal = ET.SubElement(plm, 'lines')
        crystal.set('sector', '1')
        crystal_mod = ET.SubElement(crystal, 'module')
        crystal_mod.set('sizex', 'CrystalModule_sx')
        crystal_mod.set('sizey', 'CrystalModule_sy')
        crystal_mod.set('sizez', 'CrystalModule_sz')
        crystal_mod.set('material', 'PbWO4')
        crystal_mod.set('vis', 'AnlTeal')
        crystal_wrap = ET.SubElement(crystal, 'wrapper')
        crystal_wrap.set('thickness', 'CrystalModule_wrap')
        crystal_wrap.set('material', 'Epoxy')
        crystal_wrap.set('vis', 'WhiteVis')
        crystal_pos = ET.SubElement(crystal, 'position')
        crystal_pos.set('z', 'CrystalModule_dz')
        # crystal placements (for individuals)
        # for m, (x, y) in enumerate(individual_placement(CRYSTAL_ALIGNMENT)):
        #     module = ET.SubElement(crystal, 'placement')
        #     module.set('x', '{:.3f}*mm'.format(x))
        #     module.set('y', '{:.3f}*mm'.format(y))
        #     module.set('z', '0')
        #     module.set('id', '{:d}'.format(m))
    
        # crystal placements (for lines)
        crystal.set('mirrorx', 'true')
        crystal.set('mirrory', 'true')
        for row, (begin, nmods) in enumerate(CRYSTAL_ALIGNMENT):
            line = ET.SubElement(crystal, 'line')
            line.set('x', '(CrystalModule_sx + CrystalModule_wrap)/2.')
            line.set('y', '(CrystalModule_sy + CrystalModule_wrap)*{:d}/2.'.format(row*2 + 1))
            line.set('begin', '{:d}'.format(begin - 1))
            line.set('nmods', '{:d}'.format(nmods))
    
    
        # glass
        # glass = ET.SubElement(plm, 'individuals')
        glass = ET.SubElement(plm, 'lines')
        glass.set('sector', '2')
        glass_mod = ET.SubElement(glass, 'module')
        glass_mod.set('sizex', 'GlassModule_sx')
        glass_mod.set('sizey', 'GlassModule_sy')
        glass_mod.set('sizez', 'GlassModule_sz')
        # TODO: change glass material
        glass_mod.set('material', 'PbGlass')
        glass_mod.set('vis', 'AnlBlue')
        glass_wrap = ET.SubElement(glass, 'wrapper')
        glass_wrap.set('thickness', 'GlassModule_wrap')
        glass_wrap.set('material', 'Epoxy')
        glass_wrap.set('vis', 'WhiteVis')
        glass_pos = ET.SubElement(glass, 'position')
        glass_pos.set('z', 'GlassModule_dz')
        # crystal placements (for individuals)
        # for m, (x, y) in enumerate(individual_placement(GLASS_ALIGNMENT, 41.0, 41.0)):
        #     module = ET.SubElement(glass, 'placement')
        #     module.set('x', '{:.3f}*mm'.format(x))
        #     module.set('y', '{:.3f}*mm'.format(y))
        #     module.set('z', '0')
        #     module.set('id', '{:d}'.format(m))
    
        # crystal placements (for lines)
        glass.set('mirrorx', 'true')
        glass.set('mirrory', 'true')
        for row, (begin, nmods) in enumerate(GLASS_ALIGNMENT):
            line = ET.SubElement(glass, 'line')
            line.set('x', '(GlassModule_sx + GlassModule_wrap)/2.')
            line.set('y', '(GlassModule_sy + GlassModule_wrap)*{:d}/2.'.format(row*2 + 1))
            line.set('begin', '{:d}'.format(begin - 1))
            line.set('nmods', '{:d}'.format(nmods))
    
        # readout
        readouts = ET.SubElement(data, 'readouts')
        cmt = ET.SubElement(readouts, 'comment')
        cmt.text = 'Effectively no segmentation, the segmentation is used to provide cell dimension info'
        readout = ET.SubElement(readouts, 'readout')
        readout.set('name', readout_name)
        seg = ET.SubElement(readout, 'segmentation')
        # need segmentation to provide cell dimension info
        # seg.set('type', 'NoSegmentation')
        seg.set('type', 'MultiSegmentation')
        seg.set('key', 'sector')
        crystal_seg = ET.SubElement(seg, 'segmentation')
        crystal_seg.set('name', 'CrystalSeg')
        crystal_seg.set('key_value', '1')
        crystal_seg.set('type', 'CartesianGridXY')
        crystal_seg.set('grid_size_x', 'CrystalModule_sx + CrystalModule_wrap')
        crystal_seg.set('grid_size_y', 'CrystalModule_sy + CrystalModule_wrap')
        glass_seg = ET.SubElement(seg, 'segmentation')
        glass_seg.set('name', 'GlassSeg')
        glass_seg.set('key_value', '2')
        glass_seg.set('type', 'CartesianGridXY')
        glass_seg.set('grid_size_x', 'GlassModule_sx + GlassModule_wrap')
        glass_seg.set('grid_size_y', 'GlassModule_sy + GlassModule_wrap')
        rid = ET.SubElement(readout, 'id')
        rid.text = 'system:8,sector:4,module:20,x:32:-16,y:-16'
    
    
        text = ET.tostring(data, pretty_print=True)
        with open(args.save, 'wb') as f:
            f.write(text)