Analysis/MaterialScan/README.md

+# Howto scan the material
+
+## Introduction
+Use DD4hep tool [materialScan](https://dd4hep.web.cern.ch/dd4hep/reference/materialScan_8cpp_source.html) to get the material passing through for a geatino partical (i.e. a straight line) at a given angle.
+
+## Material Scan
+
+Select the dd4hep geometry files you wish to scan, and run the given script.
+
+```bash
+# Scan the material for a specific detector
+root -l 'src/ScanMaterial.cpp("input_file.xml", "output_file.txt", bins, world_size)'
+
+# example: scan the material for the BeamPipe detector
+root -l 'src/ScanMaterial.cpp("Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01-onlyBeamPipe.xml", "BeamPipe.txt", 100, 10000)'
+```
+
+The "FILE* outFile" truncates the "stdout" output stream to the output_file.txt ("MaterialScanLogs.txt" by default), so you can't see any output. Be patient and wait for the file to be created, it may take 1-2 minutes depending on the detectors & the number of bins.
+
+*bins* is the number of bins for theta and phi, *world_size* is the size of the scanning area (mm).
+
+The output_file.txt contains all the material passing through for a straight line at the given angles with the size of bins of theta and phi.
+
+## Draw the material budget for single xml file
+
+We give an example script [options/extract_x0.py](options/extract_x0.py) to draw the material budget from the scan file generated above (e.x. MaterialScanLogs.txt).
+
+```bash
+# Draw the material budget for single xml file
+# --input_file: the scan file generated above
+# --output_file: the output file name (default: material_budget.png)
+# --detector: the detector name to be used as the title of the plot (default: all_world)
+# --bins: the bins of theta and phi (default: 100), the theta-bins is set to "bins"+1, the phi-bins is set to "bins"*2+1
+# --bias: the bias used to cut the theta range to [bias, bins-bias] (default: 0)
+python options/extract_x0.py \
+    --input_file MaterialScanLogs.txt \
+    --output_file material_budget.png \
+    --detector all_world \
+    --bins 100 \
+    --bias 0
+```
+
+Make sure the property "bins" is the same as the one used in the [Material Scan](#material-scan).
+
+The output is a png file shows the material budget at different theta and phi, and the average X0 for theta/phi.
+
+## Draw the material budget for multiple xml files
+
+For multiple xml files, we can use the script [options/extract_x0_multi.py](options/extract_x0_multi.py) to draw the material budget.
+
+```bash
+# Draw the material budget for multiple xml files
+# --input_files: the scan files generated above, files names separated with spaces
+# --labels: the labels name to be used as the legend of the plot
+# --output_file: the output file name (default: accumulated_material_budget.png)
+# --bins: the bins of theta and phi (default: 100), the theta-bins is set to "bins"+1, the phi-bins is set to "bins"*2+1
+# --bias: the bias used to cut the theta range to [bias, bins-bias] (default: 0)
+python options/extract_x0_multi.py \
+    --input_file BeamPipe.txt Lumical.txt VXD.txt FTD.txt SIT.txt \
+    --labels BeamPipe Lumical VXD FTD SIT \
+    --output_file accumulated_material_budget.png \
+    --bins 100 \
+    --bias 0
+```
+
+The output is a png file shows the average X0 for theta/phi with different labels.
\ No newline at end of file

Analysis/MaterialScan/options/extract_multi_x0.py

+import numpy as np
+import matplotlib.pyplot as plt
+import os
+import argparse
+
+def extract_x0_values(filename):
+    """Extract x0 values from a single file"""
+    x0_values = []
+    counter = 0
+    with open(filename, 'r') as file:
+        for line in file:
+            if '|     0 Average Material' in line:
+                counter += 1
+                items = [item.strip() for item in line.split('|')]
+                x0_value = float(items[1].split()[10])
+                x0_values.append(x0_value)
+                print(f"Processing file {filename}: ", counter, " x0: ", x0_value)
+    return x0_values
+
+def process_multiple_files(file_list, bins, bias):
+    """Process multiple files and return their x0 matrices"""
+    theta_bins = bins + 1
+    phi_bins = bins*2 + 1
+    range_theta = [bias, theta_bins-bias]
+    
+    matrices = []
+    for file in file_list:
+        x0_values = extract_x0_values(file)
+        x0_matrix = np.array(x0_values).reshape(theta_bins, phi_bins)[range_theta[0]:range_theta[1], :]
+        matrices.append(x0_matrix)
+    
+    return matrices
+
+def plot_accumulated_projections(matrices, output_file, bins, bias, labels):
+    """Plot accumulated projections"""
+    theta_bins = bins + 1
+    phi_bins = bins*2 + 1
+    range_theta = [bias, theta_bins-bias]
+    
+    # Create angle arrays
+    phi = np.linspace(-180, 180, phi_bins)
+    theta = np.linspace(range_theta[0]*180/theta_bins, range_theta[1]*180/theta_bins, range_theta[1]-range_theta[0])
+    
+    # Create figure
+    fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(20, 8))
+    
+    # Accumulation data
+    accumulated_theta = np.zeros_like(theta)
+    accumulated_phi = np.zeros_like(phi)
+    
+    colors = plt.cm.viridis(np.linspace(0, 1, len(matrices)))
+    
+    # Plot phi projection
+    for i, matrix in enumerate(matrices):
+        phi_projection = np.sum(matrix, axis=0) / (range_theta[1]-range_theta[0])
+        accumulated_phi += phi_projection
+        ax1.fill_between(phi, accumulated_phi, accumulated_phi - phi_projection, 
+                        label=labels[i], color=colors[i], alpha=0.6)
+    
+    # Plot theta projection
+    for i, matrix in enumerate(matrices):
+        theta_projection = np.sum(matrix, axis=1) / phi_bins
+        accumulated_theta += theta_projection
+        ax2.fill_between(theta, accumulated_theta, accumulated_theta - theta_projection,
+                        label=labels[i], color=colors[i], alpha=0.6)
+    
+    # Set phi projection plot
+    ax1.set_title('Accumulated Projection along Phi', fontsize=14, pad=10)
+    ax1.set_xlabel('Phi (degree)', fontsize=12)
+    ax1.set_ylabel('Accumulated X/X0', fontsize=12)
+    ax1.grid(True, linestyle='--', alpha=0.7)
+    ax1.legend(fontsize=10)
+    ax1.spines['top'].set_visible(False)
+    ax1.spines['right'].set_visible(False)
+    
+    # Set theta projection plot
+    ax2.set_title('Accumulated Projection along Theta', fontsize=14, pad=10)
+    ax2.set_xlabel('Theta (degree)', fontsize=12)
+    ax2.set_ylabel('Accumulated X/X0', fontsize=12)
+    ax2.grid(True, linestyle='--', alpha=0.7)
+    ax2.legend(fontsize=10)
+    ax2.spines['top'].set_visible(False)
+    ax2.spines['right'].set_visible(False)
+    
+    plt.suptitle('Detector Material Budget Accumulation Analysis', fontsize=16, y=0.95)
+    plt.tight_layout()
+    
+    plt.savefig(output_file, dpi=600, bbox_inches='tight')
+    plt.show()
+
+def main():
+    parser = argparse.ArgumentParser(description='Process multiple material scan data and generate accumulated distribution plots.')
+    parser.add_argument('--input_files', nargs='+', required=True,
+                        help='List of input files in order')
+    parser.add_argument('--labels', nargs='+', required=True,
+                        help='Labels for each input file')
+    parser.add_argument('--output_file', default='accumulated_material_budget.png',
+                        help='Output PNG file path')
+    parser.add_argument('--bins', type=int, default=100,
+                        help='theta bins is "bins"+1, phi bins is "bins"*2+1, (default: 100)')
+    parser.add_argument('--bias', type=int, default=0,
+                        help='Bias value for theta range (default: 0)')
+    
+    args = parser.parse_args()
+    
+    if len(args.input_files) != len(args.labels):
+        raise ValueError("Number of input files must match number of labels!")
+    
+    # Process all input files
+    matrices = process_multiple_files(args.input_files, args.bins, args.bias)
+    
+    # Plot accumulated graphs
+    plot_accumulated_projections(matrices, args.output_file, args.bins, args.bias, args.labels)
+
+if __name__ == '__main__':
+    main()
\ No newline at end of file

Analysis/MaterialScan/options/extract_x0.py

+import numpy as np
+import matplotlib.pyplot as plt
+import os
+import argparse
+
+def extract_x0_values(filename):
+    x0_values = []
+    counter = 0
+    with open(filename, 'r') as file:
+        for line in file:
+            if '|     0 Average Material' in line:
+                counter += 1
+                items = [item.strip() for item in line.split('|')]
+                x0_value = float(items[1].split()[10])
+                x0_values.append(x0_value)
+                print("processing: ", counter, " x0: ", x0_value)
+    return x0_values
+
+def plot_all_in_one(x0_values, output_file, bins, bias, detector):
+    """Plot all charts in one figure"""
+    theta_bins = bins + 1
+    phi_bins = bins*2 + 1
+    range_theta = [bias, theta_bins-bias]
+    # Adjust figure size and ratio, increase height to accommodate title
+    fig = plt.figure(figsize=(28, 8))
+    
+    # Adjust spacing between subplots and top margin
+    gs = fig.add_gridspec(1, 3, width_ratios=[1.2, 1, 1], wspace=0.25, top=0.85)
+    
+    # Add three subplots
+    ax1 = fig.add_subplot(gs[0])  # Heat map
+    ax2 = fig.add_subplot(gs[1])  # Phi projection
+    ax3 = fig.add_subplot(gs[2])  # Theta projection
+    
+    # Draw heat map
+    X0_matrix = np.array(x0_values).reshape(theta_bins, phi_bins)[range_theta[0]:range_theta[1], :]
+    phi = np.linspace(-180, 180, phi_bins)
+    theta = np.linspace(range_theta[0]*180/theta_bins, range_theta[1]*180/theta_bins, range_theta[1]-range_theta[0])
+    THETA, PHI = np.meshgrid(theta, phi)
+    im = ax1.pcolormesh(THETA, PHI, X0_matrix.T, shading='auto', cmap='viridis')
+    cbar = plt.colorbar(im, ax=ax1)
+    cbar.set_label('X/X0', fontsize=12)
+    ax1.set_title('Material Budget Distribution', fontsize=12, pad=10)
+    ax1.set_xlabel('Theta (angle)', fontsize=10)
+    ax1.set_ylabel('Phi (angle)', fontsize=10)
+    
+    # Draw phi projection
+    phi_projection = np.sum(X0_matrix, axis=0) / (range_theta[1]-range_theta[0])
+    ax2.plot(phi, phi_projection, 'b-', linewidth=2)
+    ax2.fill_between(phi, phi_projection, alpha=0.3)
+    ax2.set_title('Projection along Phi', fontsize=12, pad=10)
+    ax2.set_xlabel('Phi (degree)', fontsize=10)
+    ax2.set_ylabel('Average X/X0', fontsize=10)
+    ax2.grid(True, linestyle='--', alpha=0.7)
+    ax2.spines['top'].set_visible(False)
+    ax2.spines['right'].set_visible(False)
+    
+    # Draw theta projection
+    theta_projection = np.sum(X0_matrix, axis=1) / phi_bins
+    ax3.plot(theta, theta_projection, 'r-', linewidth=2)
+    ax3.fill_between(theta, theta_projection, alpha=0.3)
+    ax3.set_title('Projection along Theta', fontsize=12, pad=10)
+    ax3.set_xlabel('Theta (degree)', fontsize=10)
+    ax3.set_ylabel('Average X/X0', fontsize=10)
+    ax3.grid(True, linestyle='--', alpha=0.7)
+    ax3.spines['top'].set_visible(False)
+    ax3.spines['right'].set_visible(False)
+    
+    # Adjust main title position
+    plt.suptitle('Material Budget Analysis for {}'.format(detector), 
+                fontsize=16,  # Increase font size
+                y=0.95)  # Adjust title vertical position
+    
+    plt.tight_layout()
+    plt.savefig(output_file, dpi=600, bbox_inches='tight')
+    plt.show()
+
+def main():
+    parser = argparse.ArgumentParser(description='Process material scan data and generate plots.')
+    parser.add_argument('--input_file', help='Path to the input text file')
+    parser.add_argument('--output_file', default='material_budget.png', help='Path to the output png file')
+    parser.add_argument('--detector', default='all_world', help='Detector name')
+    parser.add_argument('--bins', type=int, default=100, help='theta bins is "bins"+1, phi bins is "bins"*2+1, (default: 100)')
+    parser.add_argument('--bias', type=int, default=0, help='Bias value for theta range (default: 0)')
+    args = parser.parse_args()
+
+    # all_world  Coil  Ecal  FTD  Hcal  Lumical  Muon  OTK  ParaffinEndcap  SIT  TPC  VXD
+    # bias detectors: BeamPipe onlyTracker
+    # Extract detector name from filename
+    x0_values = extract_x0_values(args.input_file)
+    plot_all_in_one(x0_values, args.output_file, args.bins, args.bias, args.detector)
+
+if __name__ == '__main__':
+    main() 
\ No newline at end of file

Analysis/MaterialScan/src/ScanMaterial.cpp

+#include <DD4hep/Detector.h>
+#include <DDRec/MaterialScan.h>
+#include <fstream>
+
+void ScanMaterial(string infile, string outfile="MaterialScanLogs.txt", int bins=100, double world_size=1000){
+    using namespace dd4hep;
+    using namespace dd4hep::rec;
+
+    if(infile.empty()){
+        cout << "Usage: root -l 'ScanMaterial.cpp(\"input_file.xml\", \"output_file.txt\", bins, world_size)' " << endl;
+        return;
+    }
+    
+    Detector& description = Detector::getInstance();
+    description.fromXML(infile.c_str());
+    MaterialScan scan(description);
+
+    std::ofstream clearFile(outfile.c_str(), std::ios::out | std::ios::trunc);
+    clearFile.close();
+    FILE* outFile = freopen(outfile.c_str(), "w", stdout);
+    for(int thetabin=0; thetabin<=bins; thetabin++){
+        double theta = thetabin * M_PI / bins;
+        double z = world_size*cos(theta);
+        double tranverse = world_size*sin(theta);
+        for(int phibin=-bins; phibin<=bins; phibin++){
+            double phi = phibin * M_PI / bins;
+            double x = tranverse*cos(phi);
+            double y = tranverse*sin(phi);
+            scan.print(0,0,0,x,y,z);
+        }
+    }
+    fclose(outFile);
+    return;
+}
\ No newline at end of file

Detector/DetCRD/CMakeLists.txt

@@ -8,6 +8,7 @@ gaudi_add_module(DetCRD
                          src/Calorimeter/LongCrystalBarBarrelCalorimeter32Polygon_v01.cpp
                          src/Calorimeter/LongCrystalBarEndcapCalorimeter_v01.cpp
                          src/Calorimeter/LongCrystalBarEndcapCalorimeter_v02.cpp
+                         src/Calorimeter/LongCrystalBarEndcapCalorimeter_v03.cpp
                          src/Calorimeter/CRDEcal_Short_v02.cpp
                          src/Calorimeter/CRDEcal_Endcap_Short_v01.cpp
                          src/Calorimeter/RotatedPolyhedraBarrelCalorimeter_v01_geo.cpp
@@ -15,7 +16,7 @@ gaudi_add_module(DetCRD
                          src/Calorimeter/Lumical_v01_geo.cpp
                          src/Other/Lumical_v01_geo_beampipe.cpp
                          src/Other/CRDBeamPipe_v01_geo.cpp
-                         src/Muon/Muon_Barrel_v01_03.cpp
+                         src/Muon/Muon_Barrel_v01_04.cpp
                          src/Muon/Muon_Endcap_v01_02.cpp
                          src/Tracker/SiTrackerSkewRing_v01_geo.cpp
                          src/Tracker/SiTrackerStitching_v01_geo.cpp

Detector/DetCRD/compact/CRD_common_v01/Ecal_Crystal_Endcap_v01_03.xml

+<lccdd>
+  <define>
+    <!-- <constant name="Ecal_endcap_inner_radius" value="2.5*mm"/>
+    <constant name="Ecal_endcap_inner_radius" value="2.5*mm"/> -->
+    <constant name="Ecal_endcap_nlayers" value="28"/> 
+    <constant name="Ecal_scintillator_thickness" value="10*mm"/>
+    <constant name="Ecal_deadarea_thickness" value="8.5*mm"/>
+    <constant name="Ecal_module_safety" value="0.5*mm"/>
+  </define>
+
+  <detectors>
+    <detector id="DetID_ECAL_ENDCAP" name="EcalEndcap" type="LongCrystalBarEndcapCalorimeter_v03" readout="EcalEndcapsCollection"  vis="CyanVis" calorimeterType="EMC_ENDCAP">
+      <comment>Electromagnetic Calorimeter Endcap</comment>
+
+      <envelope vis="SeeThrough">
+        <shape type="BooleanShape" operation="Subtraction" material="Air">
+          <shape type="BooleanShape" operation="Subtraction" material="Air">
+            <!-- <shape type="Tube" rmin="0.0" rmax="Ecal_endcap_outer_radius - env_safety" dz="Ecal_endcap_zmax"/>  -->
+            <!--there is a thin plane in envolop -->
+            <!-- <shape type="Tube" rmin="0.0" rmax="Ecal_endcap_outer_radius + env_safety" dz="Ecal_endcap_zmin"/> -->
+            <shape type="PolyhedraRegular" numsides="256" rmin="0.0" rmax="Hcal_barrel_inner_radius - 5*mm" dz="Ecal_endcap_zmax*2"/>
+            <shape type="PolyhedraRegular" numsides="256" rmin="0.0" rmax="Hcal_barrel_inner_radius - 5*mm" dz="Ecal_endcap_zmin*2"/>
+          </shape>
+        <!-- <shape type="Box" dx="Ecal_endcap_inner_radius" dy="Ecal_endcap_inner_radius" dz="Ecal_endcap_zmax + env_safety"/>  -->
+        <shape type="Tube" rmin="0.0" rmax="Ecal_endcap_inner_radius" dz="Ecal_endcap_zmax + env_safety"/>
+        </shape>
+        <rotation x="0" y="0" z="0"/>
+        <!-- <rotation x="0" y="0" z="11.25"/> -->
+      </envelope>
+
+      <type_flags type=" DetType_CALORIMETER + DetType_ENDCAP + DetType_EMC " />
+
+      <material name="CarbonFiber"/> 
+ 
+      <dimensions numsides="Ecal_x_module" > <!-- 0:cube 1:isosceles trapezoid 2:right trapezoid -->
+        <dimensions id="0"  module_type="0" module_number="3" x_offset="768*mm" y_offset="768*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="326*mm" dim_y="326*mm" dim_z="300*mm"/>
+        <dimensions id="1"  module_type="0" module_number="2" x_offset="1069*mm" y_offset="1069*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="276*mm" dim_y="276*mm" dim_z="300*mm"/>
+        <dimensions id="2"  module_type="0" module_number="1" x_offset="1305.5*mm" y_offset="1305.5*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="197*mm" dim_y="197*mm" dim_z="300*mm"/>
+        <dimensions id="0"  module_type="20" module_number="3" x_offset="768*mm" y_offset="768*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="326*mm" dim_y="326*mm" dim_z="300*mm"/>
+        <dimensions id="1"  module_type="20" module_number="2" x_offset="1069*mm" y_offset="1069*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="276*mm" dim_y="276*mm" dim_z="300*mm"/>
+        <dimensions id="2"  module_type="20" module_number="1" x_offset="1305.5*mm" y_offset="1305.5*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="197*mm" dim_y="197*mm" dim_z="300*mm"/>
+
+        <dimensions id="3"  module_type="1" module_number="4" x_offset="570.5*mm" y_offset="0*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="574*mm" dim_x2="495*mm" dim_y1="441*mm" dim_y2="441*mm" dim_z="300*mm" />
+        <dimensions id="4"  module_type="2" module_number="4" x_offset="570.5*mm" y_offset="0*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="574*mm" dim_x2="495*mm" dim_y1="441*mm" dim_y2="441*mm" dim_z="300*mm" />
+        <dimensions id="3"  module_type="21" module_number="4" x_offset="570.5*mm" y_offset="0*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="574*mm" dim_x2="495*mm" dim_y1="441*mm" dim_y2="441*mm" dim_z="300*mm" />
+        <dimensions id="4"  module_type="22" module_number="4" x_offset="570.5*mm" y_offset="0*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="574*mm" dim_x2="495*mm" dim_y1="441*mm" dim_y2="441*mm" dim_z="300*mm" />
+
+        <dimensions id="5"  module_type="3" module_number="4" x_offset="426.25*mm" y_offset="436.125*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="318*mm" dim_x2="357.5*mm" dim_y1="357.5*mm" dim_y2="357.5*mm" dim_z="300*mm"/>
+        <dimensions id="6"  module_type="4" module_number="4" x_offset="426.25*mm" y_offset="436.125*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="318*mm" dim_x2="357.5*mm" dim_y1="357.5*mm" dim_y2="357.5*mm" dim_z="300*mm"/>
+        <dimensions id="7"  module_type="5" module_number="4" x_offset="426.25*mm" y_offset="436.125*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="318*mm" dim_x2="357.5*mm" dim_y1="357.5*mm" dim_y2="357.5*mm" dim_z="300*mm"/>
+        <dimensions id="8"  module_type="6" module_number="4" x_offset="426.25*mm" y_offset="436.125*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="318*mm" dim_x2="357.5*mm" dim_y1="357.5*mm" dim_y2="357.5*mm" dim_z="300*mm"/>
+        <dimensions id="9"  module_type="7" module_number="4" x_offset="426.25*mm" y_offset="436.125*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="318*mm" dim_x2="357.5*mm" dim_y1="357.5*mm" dim_y2="357.5*mm" dim_z="300*mm"/>
+        <dimensions id="10"  module_type="8" module_number="4" x_offset="426.25*mm" y_offset="436.125*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="318*mm" dim_x2="357.5*mm" dim_y1="357.5*mm" dim_y2="357.5*mm" dim_z="300*mm"/>
+        <dimensions id="11"  module_type="9" module_number="4" x_offset="426.25*mm" y_offset="436.125*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="318*mm" dim_x2="357.5*mm" dim_y1="357.5*mm" dim_y2="357.5*mm" dim_z="300*mm"/>
+        <dimensions id="12"  module_type="10" module_number="4" x_offset="426.25*mm" y_offset="436.125*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x1="318*mm" dim_x2="357.5*mm" dim_y1="357.5*mm" dim_y2="357.5*mm" dim_z="300*mm"/>
+
+        <dimensions id="13"  module_type="11" module_number="1" x_offset="426.25*mm" y_offset="426.25*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="178.75*mm" dim_y="139.25*mm" dim_z="300*mm"/>
+        <dimensions id="14"  module_type="12" module_number="1" x_offset="426.25*mm" y_offset="426.25*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="178.75*mm" dim_y="139.25*mm" dim_z="300*mm"/>
+        <dimensions id="15"  module_type="13" module_number="1" x_offset="426.25*mm" y_offset="426.25*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="178.75*mm" dim_y="139.25*mm" dim_z="300*mm"/>
+        <dimensions id="16"  module_type="14" module_number="1" x_offset="426.25*mm" y_offset="426.25*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="178.75*mm" dim_y="139.25*mm" dim_z="300*mm"/>
+        <dimensions id="17"  module_type="15" module_number="1" x_offset="426.25*mm" y_offset="426.25*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="178.75*mm" dim_y="139.25*mm" dim_z="300*mm"/>
+        <dimensions id="18"  module_type="16" module_number="1" x_offset="426.25*mm" y_offset="426.25*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="178.75*mm" dim_y="139.25*mm" dim_z="300*mm"/>
+        <dimensions id="19"  module_type="17" module_number="1" x_offset="426.25*mm" y_offset="426.25*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="178.75*mm" dim_y="139.25*mm" dim_z="300*mm"/>
+        <dimensions id="20"  module_type="18" module_number="1" x_offset="426.25*mm" y_offset="426.25*mm" z_offset="Ecal_endcap_zmin+Ecal_endcap_thickness/2" dim_x="178.75*mm" dim_y="139.25*mm" dim_z="300*mm"/>
+        
+      </dimensions>
+
+      <layer repeat="Ecal_endcap_nlayers" vis="CyanVis" material="G4_BGO" thickness = "Ecal_scintillator_thickness">
+        <slice material="G4_BGO" thickness = "Ecal_scintillator_thickness" sensitive = "yes"   limits="cal_limits"  vis="CyanVis"   />
+      </layer>
+    </detector>
+  </detectors>
+
+  <readouts>
+    <readout name="EcalEndcapsCollection">
+      <segmentation type="NoSegmentation"/>
+      <id>system:5,module:1,part:7,stave:7,type:4,dlayer:4,slayer:1,bar:7</id>
+    </readout>
+  </readouts>
+
+</lccdd>
+
+
+
+

Detector/DetCRD/compact/CRD_common_v01/Muon_Barrel_v01_04.xml

+<?xml version="1.0" encoding="UTF-8"?>
+<lccdd>
+ 
+  <info name="Muon_Barrel"
+	title="Test with A Single Muon Barrel"
+	author="Zibing Bai"
+	url="http://cepcgit.ihep.ac.cn"
+	status="development"
+	version="v01">
+    <comment>Test with A Single Muon Barrel</comment>        
+  </info>
+
+  <define>
+    <!--Muon Barrel-->
+    <constant name="Muon_barrel_superlayer_num" value="6"/>
+    <constant name="Muon_barrel_strip_num_0" value="26"/>
+    <constant name="Muon_barrel_strip_num_1" value="38"/>
+    <constant name="Muon_barrel_strip_num_2" value="50"/>
+    <constant name="Muon_barrel_strip_num_3" value="62"/>
+    <constant name="Muon_barrel_strip_num_4" value="74"/>
+    <constant name="Muon_barrel_strip_num_5" value="86"/>
+    <!--constant name="Muon_barrel_strip_num_6" value="89"/>
+    <constant name="Muon_barrel_strip_num_7" value="101"/-->
+    <constant name="Muon_barrel_strip_num_fixed_0" value="106"/>
+    <constant name="Muon_barrel_strip_num_fixed_1" value="115"/>
+
+    <!--constant name="Muon_barrel_iron_x1" value="Muon_standard_scale"/-->
+    <constant name="Muon_barrel_iron_y" value="Muon_total_length"/>
+    <constant name="Muon_barrel_iron_z" value="Muon_standard_scale"/>
+    <constant name="Muon_barrel_iron_posx" value="-1*Muon_standard_scale"/>
+
+    <constant name="Muon_barrel_barrel_y" value="Muon_barrel_iron_y"/>
+    <constant name="Muon_barrel_barrel_posy" value="0.5*Muon_barrel_barrel_y"/>
+
+    <constant name="Muon_barrel_superlayer_init" value="-35*cm"/>
+    <constant name="Muon_barrel_superlayer_gap" value="14*cm"/>
+    <constant name="Muon_barrel_superlayer_endcap_gap" value="10*cm"/>
+    <constant name="Muon_barrel_superlayer_air_gap" value="1*cm"/>
+
+    <constant name="Muon_barrel_superlayer_length_0" value="4275*mm"/>
+    <constant name="Muon_barrel_superlayer_length_1" value="4625*mm"/>
+
+    <constant name="Muon_barrel_strip_width_fixed" value="20*mm"/>
+
+    <constant name="Muon_barrel_superlayer_y" value="2*Muon_strip_y+Muon_barrel_superlayer_air_gap"/>
+    <!--constant name="Muon_barrel_superlayer_z" value="Muon_strip_z+2*Muon_strip_surf+2*Muon_barrel_superlayer_air_gap"/-->
+
+    <!--Checkout-->
+    <!--constant name="Muon_barrel_inner_radius" value="4245*mm"/-->
+    <constant name="Muon_barrel_barrel_num" value="2"/>
+    <constant name="Muon_barrel_iron_part_num" value="12"/>
+
+  </define>
+
+  <detectors>
+    <detector id="DetID_MUON" name="MuonBarrel" type="Muon_Barrel_v01_04" readout="MuonBarrelCollection" vis="WhiteVis">
+      <position x="0" y="0" z="0"/>
+      <barrel id="Muon_barrel_iron_part_num" name="Muon_barrel_barrel" type="Muon_barrel_barrel" vis="SeeThrough">
+        <position x="0" y="Muon_barrel_barrel_posy" z="0"/>
+        <iron id="0" name="Muon_barrel_iron_part" type="Muon_barrel_iron_part" vis="GrayVis" material="Iron">
+          <material name="Iron"/>
+          <position x="Muon_barrel_iron_posx" y="0"/>
+          <dimensions x1="0.5*Muon_barrel_iron_x1" y1="0.5*Muon_barrel_iron_y" y2="0.5*Muon_barrel_iron_y" dz="0.5*Muon_barrel_iron_z"/>
+          <superlayer id="Muon_barrel_superlayer_num" name="Muon_barrel_superlayer" type="Muon_barrel_superlayer" vis="BlueVis" material="Air">
+              <stripe id="0" name="Muon_stripe" type="Muon_stripe" vis="GreenVis" material="BC420">
+                <material name="BC420"/>
+                <dimensions dx="0.5*Muon_strip_x"  dy="0.5*Muon_strip_y"  dz="0.5*Muon_strip_z+Muon_strip_SiPM_z"/>
+              </stripe>
+          </superlayer>
+        </iron>
+      </barrel>
+    </detector>
+  </detectors>
+
+  <readouts>
+    <readout name="MuonBarrelCollection">
+      <segmentation type="NoSegmentation"/>
+      <id>system:5,Env:5,Fe:5,Superlayer:5,Layer:5,Stripe:9</id>
+    </readout>
+  </readouts>
+
+</lccdd>

Detector/DetCRD/compact/TDR_CaloTest/TDR_longEcalEndcap.xml

+<?xml version="1.0" encoding="UTF-8"?>
+<lccdd xmlns:compact="http://www.lcsim.org/schemas/compact/1.0"
+       xmlns:xs="http://www.w3.org/2001/XMLSchema"
+       xs:noNamespaceSchemaLocation="http://www.lcsim.org/schemas/compact/1.0/compact.xsd">
+  <info name="TDR_o1_v01"
+        title="CepC reference detctor for TDR"
+        author=""
+        url="http://cepc.ihep.ac.cn"
+        status="developing"
+        version="v01">
+    <comment>CepC reference detector simulation models used for TDR </comment>
+  </info>
+  
+  <includes>
+    <gdmlFile  ref="${DD4hepINSTALL}/DDDetectors/compact/elements.xml"/>
+    <gdmlFile  ref="../CRD_common_v02/materials.xml"/>
+  </includes>
+  
+  <define>
+    <constant name="world_size" value="10*m"/>
+    <constant name="world_x" value="world_size"/>
+    <constant name="world_y" value="world_size"/>
+    <constant name="world_z" value="world_size"/>
+
+    <include ref="${DD4hepINSTALL}/DDDetectors/compact/detector_types.xml"/>
+  </define>
+
+  <include ref="./TDR_Dimensions_v01_01.xml"/>
+
+  <!--TODO: vertex cooling-->
+  <!-- <include ref="../CRD_common_v02/Beampipe_v01_03.xml"/> -->
+  <!--preliminary vertex and tracker, to update/-->
+  <!-- <include ref="../CRD_common_v02/VXD_StaggeredLadder_v02_01.xml"/> -->
+  <!-- <include ref="../CRD_common_v02/FTD_SkewRing_v01_05.xml"/> -->
+  <!-- <include ref="../CRD_common_v02/SIT_SimplePixel_v01_03.xml"/> -->
+  <!--include ref="../CRD_common_v01/TPC_Simple_v10_02.xml"/-->
+  <!--use 10 rows clustering version/-->
+  <!-- <include ref="../CRD_common_v02/TPC_ModularEndcap_o1_v02.xml"/> -->
+  <!-- <include ref="../CRD_common_v01/SET_SimplePixel_v01_01.xml"/> -->
+
+  <!-- <include ref="../CRD_common_v01/Ecal_Crystal_Barrel_v01_02.xml"/> -->
+  <include ref="../CRD_common_v01/Ecal_Crystal_Endcap_v01_02.xml"/> 
+  <!-- <include ref="../CRD_common_v01/SHcalGlass_Barrel_v05.xml"/> -->
+  <!-- <include ref="../CRD_common_v01/SHcalGlass_Endcaps_v01.xml"/> -->
+
+  <!--Lumical to update-->
+  <!-- <include ref="../CRD_common_v01/Lumical_o1_v01.xml"/> -->
+  <!--preliminary Magnet, to update/-->
+  <!-- <include ref="../CRD_common_v02/Coil_Simple_v01_02.xml"/> -->
+  <!--preliminary Muon, obselete/-->
+  <!--include ref="../CRD_common_v02/Yoke_Polyhedra_Barrel_v01_01.xml"/>
+  <include ref="../CRD_common_v02/Yoke_Polyhedra_Endcaps_v01_01.xml"/-->
+ 
+  <!--muon detector-->
+  <!-- <include ref="../CRD_common_v01/Muon_Barrel_v01_01.xml"/> -->
+  <!-- <include ref="../CRD_common_v01/Muon_Endcap_v01_01.xml"/> -->
+ 
+  <fields>
+    <field name="InnerSolenoid" type="solenoid"
+           inner_field="Field_nominal_value"
+           outer_field="0"
+           zmax="SolenoidCoil_half_length"
+           inner_radius="SolenoidCoil_center_radius"
+           outer_radius="Solenoid_outer_radius">
+    </field>
+    <field name="OuterSolenoid" type="solenoid"
+           inner_field="0"
+           outer_field="Field_outer_nominal_value"
+           zmax="SolenoidCoil_half_length"
+           inner_radius="Solenoid_outer_radius"
+           outer_radius="Yoke_barrel_inner_radius">
+    </field>
+  </fields>
+
+</lccdd>

Detector/DetCRD/compact/TDR_CaloTest/TDR_onlyLongEcal.xml

+<?xml version="1.0" encoding="UTF-8"?>
+<lccdd xmlns:compact="http://www.lcsim.org/schemas/compact/1.0"
+       xmlns:xs="http://www.w3.org/2001/XMLSchema"
+       xs:noNamespaceSchemaLocation="http://www.lcsim.org/schemas/compact/1.0/compact.xsd">
+  <info name="TDR_o1_v01"
+        title="CepC reference detctor for TDR"
+        author=""
+        url="http://cepc.ihep.ac.cn"
+        status="developing"
+        version="v01">
+    <comment>CepC reference detector simulation models used for TDR </comment>
+  </info>
+  
+  <includes>
+    <gdmlFile  ref="${DD4hepINSTALL}/DDDetectors/compact/elements.xml"/>
+    <gdmlFile  ref="../CRD_common_v02/materials.xml"/>
+  </includes>
+  
+  <define>
+    <constant name="world_size" value="10*m"/>
+    <constant name="world_x" value="world_size"/>
+    <constant name="world_y" value="world_size"/>
+    <constant name="world_z" value="world_size"/>
+
+    <include ref="${DD4hepINSTALL}/DDDetectors/compact/detector_types.xml"/>
+  </define>
+
+  <include ref="./TDR_Dimensions_v01_01.xml"/>
+
+  <!--TODO: vertex cooling-->
+  <!--include ref="../CRD_common_v02/Beampipe_v01_03.xml"/-->
+  <!--preliminary vertex and tracker, to update/-->
+  <!--include ref="../CRD_common_v02/VXD_StaggeredLadder_v02_01.xml"/>
+  <include ref="../CRD_common_v02/FTD_SkewRing_v01_05.xml"/>
+  <include ref="../CRD_common_v02/SIT_SimplePixel_v01_03.xml"/-->
+  <!--include ref="../CRD_common_v01/TPC_Simple_v10_02.xml"/-->
+  <!--use 10 rows clustering version/-->
+  <!--include ref="../CRD_common_v02/TPC_ModularEndcap_o1_v02.xml"/>
+  <include ref="../CRD_common_v01/SET_SimplePixel_v01_01.xml"/-->
+
+  <include ref="../CRD_common_v01/Ecal_Crystal_Barrel_v01_02.xml"/>
+  <include ref="../CRD_common_v01/Ecal_Crystal_Endcap_v01_01.xml"/>
+  <!--include ref="../CRD_common_v01/SHcalSc04_Barrel_v04_01.xml"/-->
+  <!--include ref="../CRD_common_v01/SHcalSc04_Endcaps_v02.xml"/-->
+
+  <!--Lumical to update-->
+  <!--include ref="../CRD_common_v01/Lumical_o1_v01.xml"/-->
+  <!--preliminary Magnet, to update/-->
+  <!--include ref="../CRD_common_v02/Coil_Simple_v01_02.xml"/-->
+  <!--preliminary Muon, obselete/-->
+  <!--include ref="../CRD_common_v02/Yoke_Polyhedra_Barrel_v01_01.xml"/>
+  <include ref="../CRD_common_v02/Yoke_Polyhedra_Endcaps_v01_01.xml"/-->
+ 
+  <!--muon detector-->
+  <!--include ref="../CRD_common_v01/Muon_Barrel_v01_01.xml"/-->
+  <!--include ref="../CRD_common_v01/Muon_Endcap_v01_01.xml"/-->
+ 
+  <fields>
+    <field name="InnerSolenoid" type="solenoid"
+           inner_field="0"
+           outer_field="0"
+           zmax="SolenoidCoil_half_length"
+           inner_radius="SolenoidCoil_center_radius"
+           outer_radius="Solenoid_outer_radius">
+    </field>
+    <field name="OuterSolenoid" type="solenoid"
+           inner_field="0"
+           outer_field="0"
+           zmax="SolenoidCoil_half_length"
+           inner_radius="Solenoid_outer_radius"
+           outer_radius="Yoke_barrel_inner_radius">
+    </field>
+  </fields>
+
+</lccdd>

Detector/DetCRD/compact/TDR_o1_v01/TDR_Dimensions_v01_01.xml

@@ -162,7 +162,7 @@
     <constant name="Ecal_barrel_symmetry"     value="32"/>
     <constant name="Ecal_Tpc_gap"             value="Ecal_barrel_inner_radius-TPC_outer_radius"/>
 
-    <constant name="Ecal_endcap_inner_radius" value="400*mm"/>
+    <constant name="Ecal_endcap_inner_radius" value="350*mm"/>
     <constant name="Ecal_endcap_outer_radius" value="Ecal_barrel_outer_radius"/>
     <constant name="Ecal_endcap_zmin"         value="2930*mm"/>
     <constant name="Ecal_endcap_thickness"    value="Ecal_barrel_thickness"/>
@@ -263,7 +263,7 @@
     <!--standard scale-->
     <constant name="Muon_standard_scale" value="94*cm"/>    
     <!--constant name="Muon_total_length" value="8950*mm"/--><!--overlap with coil and hcal-->
-    <constant name="Muon_total_length" value="Yoke_endcap_zmin*2"/>
+    <constant name="Muon_total_length" value="9150*mm"/>
 
     <!--Muon Barrel>
     <constant name="Muon_barrel_barrel_num" value="2"/>

Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01-onlyEcal.xml

+<?xml version="1.0" encoding="UTF-8"?>
+<lccdd xmlns:compact="http://www.lcsim.org/schemas/compact/1.0"
+       xmlns:xs="http://www.w3.org/2001/XMLSchema"
+       xs:noNamespaceSchemaLocation="http://www.lcsim.org/schemas/compact/1.0/compact.xsd">
+  <info name="TDR_o1_v01"
+        title="CepC reference detctor for TDR"
+        author=""
+        url="http://cepc.ihep.ac.cn"
+        status="developing"
+        version="v01">
+    <comment>CepC reference detector simulation models used for TDR </comment>
+  </info>
+  
+  <includes>
+    <gdmlFile  ref="${DD4hepINSTALL}/DDDetectors/compact/elements.xml"/>
+    <gdmlFile  ref="../CRD_common_v02/materials.xml"/>
+  </includes>
+  
+  <define>
+    <constant name="world_size" value="10*m"/>
+    <constant name="world_x" value="world_size"/>
+    <constant name="world_y" value="world_size"/>
+    <constant name="world_z" value="world_size"/>
+
+    <include ref="${DD4hepINSTALL}/DDDetectors/compact/detector_types.xml"/>
+  </define>
+
+  <include ref="./TDR_Dimensions_v01_01.xml"/>
+
+  <include ref="../CRD_common_v01/Ecal_Crystal_Barrel_v01_02.xml"/>
+  <!--preliminary EcalEndcaps/-->
+  <include ref="../CRD_common_v02/EcalEndcaps_Polyhedra_v01_01.xml"/>
+ 
+  <fields>
+    <field name="InnerSolenoid" type="solenoid"
+           inner_field="Field_nominal_value"
+           outer_field="0"
+           zmax="SolenoidCoil_half_length"
+           inner_radius="SolenoidCoil_center_radius"
+           outer_radius="Solenoid_outer_radius">
+    </field>
+    <field name="OuterSolenoid" type="solenoid"
+           inner_field="0"
+           outer_field="Field_outer_nominal_value"
+           zmax="SolenoidCoil_half_length"
+           inner_radius="Solenoid_outer_radius"
+           outer_radius="Yoke_barrel_inner_radius">
+    </field>
+  </fields>
+
+</lccdd>

Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01.xml

@@ -42,7 +42,7 @@
   <!--include ref="../CRD_common_v01/OTKEndcap_v01_01.xml"/-->
 
   <include ref="../CRD_common_v01/Ecal_Crystal_Barrel_v01_02.xml"/>
-  <include ref="../CRD_common_v01/Ecal_Crystal_Endcap_v01_02.xml"/>
+  <include ref="../CRD_common_v01/Ecal_Crystal_Endcap_v01_03.xml"/>
   <include ref="../CRD_common_v01/SHcalGlass_Barrel_v05.xml"/>
   <include ref="../CRD_common_v01/SHcalGlass_Endcaps_v01.xml"/>
 
@@ -55,7 +55,7 @@
   <include ref="../CRD_common_v02/Yoke_Polyhedra_Endcaps_v01_01.xml"/-->
  
   <!--muon detector-->
-  <include ref="../CRD_common_v01/Muon_Barrel_v01_03.xml"/>
+  <include ref="../CRD_common_v01/Muon_Barrel_v01_04.xml"/>
   <include ref="../CRD_common_v01/Muon_Endcap_v01_02.xml"/>
   <include ref="../CRD_common_v01/ParaffinEndcap_v01_01.xml"/>
  

Detector/DetCRD/scripts/TDR_o1_v01/tracking.py

@@ -42,7 +42,7 @@ tracksystemsvc = TrackSystemSvc("TrackSystemSvc")
 from Configurables import SimplePIDSvc
 pidsvc = SimplePIDSvc("SimplePIDSvc")
 cepcswdatatop = "/cvmfs/cepcsw.ihep.ac.cn/prototype/releases/data/latest"
-pidsvc.ParFile = os.path.join(cepcswdatatop, "CEPCSWData/offline-data/Service/SimplePIDSvc/data/dNdx_TPC.root")
+pidsvc.ParFile = os.path.join(cepcswdatatop, "CEPCSWData/offline-data/Service/SimplePIDSvc/data/tdr25.1.0/dNdx_TPC.root")
 
 from Configurables import PodioInput
 podioinput = PodioInput("PodioReader", collections=[

Detector/DetCRD/src/Muon/Muon_Barrel_v01_04.cpp

+#include "DD4hep/DetFactoryHelper.h"
+#include "DD4hep/OpticalSurfaces.h"
+#include "XML/Utilities.h"
+#include "DDRec/DetectorData.h"
+#include "DDSegmentation/Segmentation.h"
+#include <cmath>
+
+#define MYDEBUG(x) std::cout << __FILE__ << ":" << __LINE__ << ": " << x << std::endl;
+#define MYDEBUGVAL(x) std::cout << __FILE__ << ":" << __LINE__ << ": " << #x << ": " << x << std::endl;
+
+using dd4hep::rec::LayeredCalorimeterData;
+static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
+                                     xml_h e,
+                                     dd4hep::SensitiveDetector sens) {
+
+    xml_det_t x_det = e;
+
+    std::string det_name = x_det.nameStr();
+    std::string det_type = x_det.typeStr();
+    MYDEBUGVAL(det_name);
+    MYDEBUGVAL(det_type);
+    xml_dim_t   pos    (x_det.child(_U(position)));
+
+    dd4hep::DetElement sdet(det_name, x_det.id());
+
+    dd4hep::Assembly envelope(det_name);
+
+    dd4hep::xml::setDetectorTypeFlag( e, sdet ) ;
+
+    dd4hep::Volume motherVol = theDetector.pickMotherVolume(sdet);
+
+    sens.setType("muonbarrel");
+
+//    dd4hep::Readout readout = sens.readout();
+//    dd4hep::Segmentation seg = readout.segmentation();
+    xml_coll_t dcEnv(x_det,Unicode("barrel"));
+    xml_comp_t x_env = dcEnv;
+    xml_coll_t dcFe(x_env,Unicode("iron"));
+    xml_comp_t x_Fe = dcFe;
+
+    std::string Fe_name = x_Fe.nameStr();
+    dd4hep::Material Fe_mat(theDetector.material(x_Fe.materialStr()));
+    xml_dim_t Fe_pos(x_Fe.child(_U(position)));
+    xml_dim_t Fe_dim(x_Fe.child(_U(dimensions)));
+    //double Fe_halfX1 = 0.5 * theDetector.constant<double>("Muon_standard_scale");
+   // double Fe_halfX2 = Fe_halfX1 + theDetector.constant<double>("Muon_barrel_iron_z") * std::sqrt(3);
+
+   // dd4hep::Trd2 Fe_solid(Fe_halfX1,Fe_halfX2,Fe_dim.y1(),Fe_dim.y2(),Fe_dim.dz());
+   // dd4hep::Volume Fe_vol(Fe_name, Fe_solid, Fe_mat);
+   // Fe_vol.setVisAttributes(theDetector.visAttributes(x_Fe.visStr()));
+
+    xml_coll_t dcSuperlayer(x_Fe,Unicode("superlayer"));
+    xml_comp_t x_superlayer = dcSuperlayer;
+
+    int strip_type_num_total = 8;
+    dd4hep::Box strip_solid[strip_type_num_total], strip_solid_fixed[strip_type_num_total];//, surface_solid[strip_type_num_total], BC420_solid[strip_type_num_total];
+    //dd4hep::Tube cut3_solid[strip_type_num_total];
+    dd4hep::Volume strip_vol[strip_type_num_total], strip_vol_fixed[strip_type_num_total];//, surface_vol[strip_type_num_total], BC420_vol[strip_type_num_total], cut3_vol[strip_type_num_total];
+
+    xml_coll_t dcStrip(x_superlayer,Unicode("stripe"));
+    xml_comp_t x_strip = dcStrip;
+    std::string strip_name = x_strip.nameStr();
+    dd4hep::Material strip_mat(theDetector.material(x_strip.materialStr()));
+    xml_dim_t strip_dim(x_strip.child(_U(dimensions)));
+    double fixed_width = theDetector.constant<double>("Muon_barrel_strip_width_fixed");
+
+    for (int i = 0; i < strip_type_num_total; i++ )
+    {
+      int type_num = i;
+      double strip_halfZ;
+      if ( i < 6 )
+      {
+        std::string num_name = "Muon_barrel_strip_num" + dd4hep::_toString(i,"_%d");
+        int strip_num = theDetector.constant<int>(num_name);
+        strip_halfZ = 0.5 * strip_num * theDetector.constant<double>("Muon_strip_x") + theDetector.constant<double>("Muon_strip_surf");
+      }
+
+      if ( i >= 6)
+      {
+        int fixed_num = i - 6;
+        std::string num_name = "Muon_barrel_strip_num_fixed" + dd4hep::_toString( fixed_num, "_%d");
+        int strip_num = theDetector.constant<int>(num_name);
+        strip_halfZ = 0.5 * fixed_width + 0.5 * strip_num * theDetector.constant<double>("Muon_strip_x") + theDetector.constant<double>("Muon_strip_surf");
+      }
+      //double surface_halfZ, BC420_halfZ, fiber_halfZ, cut_halfZ;
+      
+      //surface_halfZ = BC420_halfZ = fiber_halfZ = cut_halfZ = strip_halfZ - theDetector.constant<double>("Muon_strip_surf");   
+
+      std::string strip_name = x_strip.nameStr() + dd4hep::_toString( i, "_%d");
+      std::string strip_name_fixed = x_strip.nameStr() + dd4hep::_toString( i, "_%d");
+
+      strip_solid[type_num] = dd4hep::Box(strip_dim.dx(),strip_dim.dy(),strip_halfZ);
+      strip_solid_fixed[type_num] = dd4hep::Box(fixed_width, strip_dim.dy(), strip_halfZ);
+      strip_vol[type_num] = dd4hep::Volume(strip_name,strip_solid[type_num],strip_mat);
+      strip_vol_fixed[type_num] = dd4hep::Volume(strip_name_fixed,strip_solid_fixed[type_num],strip_mat);
+      strip_vol[type_num].setVisAttributes(theDetector.visAttributes(x_strip.visStr()));
+      //strip_vol[type_num].setSensitiveDetector(sens);
+      /*surface_solid[type_num] = dd4hep::Box(surface_dim.dx(),surface_dim.dy(),surface_halfZ);
+      BC420_solid[type_num] = dd4hep::Box(BC420_dim.dx(),BC420_dim.dy(),BC420_halfZ);
+      cut3_solid[type_num] = dd4hep::Tube(cut3_dim.rmin(),cut3_dim.rmax(),cut_halfZ);
+      surface_vol[type_num] = dd4hep::Volume(surface_name, surface_solid[type_num], surface_mat);
+      surface_vol[type_num].setVisAttributes(theDetector.visAttributes(surface_vis));
+
+      cut3_vol[type_num] = dd4hep::Volume(cut3_name, cut3_solid[type_num], cut3_mat);
+      cut3_vol[type_num].setVisAttributes(theDetector.visAttributes(cut3_vis));
+
+      BC420_vol[type_num] = dd4hep::Volume(BC420_name,BC420_solid[type_num],BC420_mat);
+      BC420_vol[type_num].setVisAttributes(theDetector.visAttributes(BC420_vis));
+      BC420_vol[type_num].setSensitiveDetector(sens); */
+    }
+
+    for(int i0 = 0; i0 < x_env.id(); i0++)
+    {
+      std::string env_name = x_env.nameStr() + dd4hep::_toString(i0,"_%d");
+      xml_dim_t env_pos(x_env.child(_U(position)));
+
+      /*std::string Fe_name = x_Fe.nameStr() + dd4hep::_toString(i0,"_%d");
+      dd4hep::Material Fe_mat(theDetector.material(x_Fe.materialStr()));
+      xml_dim_t Fe_pos(x_Fe.child(_U(position)));
+      xml_dim_t Fe_dim(x_Fe.child(_U(dimensions)));*/
+      double Fe_halfX1 = 0.5 * theDetector.constant<double>("Muon_standard_scale");
+      double Fe_halfX2 = Fe_halfX1 + theDetector.constant<double>("Muon_barrel_iron_z") * std::sqrt(3);
+      double Fe_posZ = -1 * theDetector.constant<double>("Muon_standard_scale") * ( 2.5 + 1.5 * std::sqrt(3) );
+      dd4hep::Assembly env_vol(env_name);
+      double env_rot = i0 * 360 * dd4hep::degree / x_env.id();
+      dd4hep::Transform3D env_transform(dd4hep::Rotation3D(dd4hep::RotationY(env_rot)),dd4hep::Position(env_pos.x(), 0,env_pos.z()));
+
+      dd4hep::Trd2 Fe_solid(Fe_halfX1,Fe_halfX2,Fe_dim.y1(),Fe_dim.y2(),Fe_dim.dz());
+      dd4hep::Volume Fe_vol(Fe_name, Fe_solid, Fe_mat);
+      Fe_vol.setVisAttributes(theDetector.visAttributes(x_Fe.visStr()));
+      dd4hep::Transform3D Fe_transform(dd4hep::Rotation3D(),dd4hep::Position(Fe_pos.x(),Fe_pos.y(),Fe_posZ));
+      for(int i1 = 0; i1 < theDetector.constant<int>("Muon_barrel_barrel_num"); i1++ )
+      {
+        xml_coll_t dcSuperlayer(x_Fe,Unicode("superlayer"));
+        xml_comp_t x_superlayer = dcSuperlayer;
+        for(int i2 = 0; i2 < x_superlayer.id(); i2++)  // FIXME:: 8 layers start from 0
+        {
+          std::string superlayer_name = x_superlayer.nameStr() + dd4hep::_toString(i1,"_%d") + dd4hep::_toString(i2,"_%d");
+          std::string num_name = "Muon_barrel_strip_num" + dd4hep::_toString(i2,"_%d");
+
+          int strip_num = theDetector.constant<int>(num_name);
+          int superlayer_posy_num = ( (2 * i1 - 1 ) * (2 * ( i2 % 2 ) - 1 ) + 1 ) / 2;
+          std::string superlayer_posy_name = "Muon_barrel_strip_num_fixed" + dd4hep::_toString( superlayer_posy_num, "_%d");
+          std::string superlayer_length_name = "Muon_barrel_superlayer_length" + dd4hep::_toString( superlayer_posy_num, "_%d");
+          double superlayer_length = theDetector.constant<double>(superlayer_length_name);
+
+          double superlayer_halfX = 0.5 * strip_num * theDetector.constant<double>("Muon_strip_x") + theDetector.constant<double>("Muon_barrel_superlayer_air_gap");
+          double superlayer_halfY = 0.5 * theDetector.constant<double>("Muon_barrel_superlayer_y");
+          double superlayer_halfZ = 0.5 * superlayer_length;
+          
+          double superlayer_posy = ( 2 * i1 - 1 ) * ( 0.5 * theDetector.constant<double>("Muon_total_length") - theDetector.constant<double>("Muon_barrel_superlayer_endcap_gap") - superlayer_halfZ);
+          double superlayer_posZ = i2 * theDetector.constant<double>("Muon_barrel_superlayer_gap") + theDetector.constant<double>("Muon_barrel_superlayer_init");
+          dd4hep::Assembly superlayer_vol(superlayer_name);
+          dd4hep::Transform3D superlayer_transform(dd4hep::Rotation3D(dd4hep::RotationX(90*dd4hep::degree)),dd4hep::Position(0, superlayer_posy, superlayer_posZ));
+
+          for (int i3 = 0; i3 < 2; i3++ )
+          {
+            int num;
+            int type;
+            if ( i3 == 0 )
+            {
+              num = strip_num;
+            }
+            if ( i3 == 1 )
+            {
+              num = theDetector.constant<int>(superlayer_posy_name);
+              strip_vol_fixed[i2].setSensitiveDetector(sens);
+              dd4hep::Transform3D strip_transform_fixed(dd4hep::Rotation3D(dd4hep::RotationY(90*dd4hep::degree)), dd4hep::Position(0, -0.5 * theDetector.constant<double>("Muon_strip_y"), theDetector.constant<double>("Muon_strip_x") * ( num - 0.5 ) - 0.5 * theDetector.constant<int>(superlayer_posy_name) * theDetector.constant<double>("Muon_strip_x")));
+              dd4hep::PlacedVolume strip_place_fixed = superlayer_vol.placeVolume(strip_vol_fixed[i2],strip_transform_fixed);
+              strip_place_fixed.addPhysVolID("Stripe", num+1).addPhysVolID("Layer",i3+1).addPhysVolID("Superlayer",i2+1).addPhysVolID("Env",i1+1);
+            }
+            for ( int i4 = 0; i4 < num; i4++ )
+            {
+              double strip_posX, strip_posY, strip_posZ;
+              dd4hep::Rotation3D strip_rot;
+              if ( i3 == 0 )
+              {
+                type = superlayer_posy_num + 6;
+                strip_posX = theDetector.constant<double>("Muon_strip_x") * ( i4 + 0.5 * (1 - strip_num));
+                strip_posY = 0.5 * theDetector.constant<double>("Muon_strip_y");
+                strip_posZ = 0;
+              }
+              if ( i3 == 1 )
+              {
+                type = i2;
+                strip_rot = dd4hep::Rotation3D(dd4hep::RotationY( 90 * dd4hep::degree ));
+                strip_posX = 0;
+                strip_posY = -0.5 * theDetector.constant<double>("Muon_strip_y");
+                strip_posZ = theDetector.constant<double>("Muon_strip_x") * ( i4 + 0.5 ) - 0.5 * theDetector.constant<int>(superlayer_posy_name) * theDetector.constant<double>("Muon_strip_x") - 0.5 * fixed_width;
+              }
+              strip_vol[type].setSensitiveDetector(sens);
+
+              dd4hep::Transform3D strip_transform(strip_rot,dd4hep::Position(strip_posX,strip_posY,strip_posZ));
+	      dd4hep::PlacedVolume strip_place = superlayer_vol.placeVolume(strip_vol[type],strip_transform);
+              strip_place.addPhysVolID("Stripe",i4+1).addPhysVolID("Layer",i3+1).addPhysVolID("Superlayer",i2+1).addPhysVolID("Env",i1+1);
+            }
+          }
+
+          dd4hep::PlacedVolume Superlayer_place = Fe_vol.placeVolume(superlayer_vol,superlayer_transform);
+        }
+      }
+      dd4hep::PlacedVolume Fe_place = env_vol.placeVolume(Fe_vol,Fe_transform);
+      Fe_place.addPhysVolID("Fe",i0+1);
+      envelope.placeVolume(env_vol,env_transform);
+    }
+    dd4hep::Transform3D pv(dd4hep::Rotation3D(dd4hep::RotationX(90*dd4hep::degree)),dd4hep::Position(0,0,0));
+    dd4hep::PlacedVolume phv = motherVol.placeVolume(envelope,pv);
+    phv.addPhysVolID("system",x_det.id());
+    sdet.setPlacement(phv);
+
+    MYDEBUG("create_detector DONE. ");
+    return sdet;
+}
+
+DECLARE_DETELEMENT(Muon_Barrel_v01_04, create_detector)

Detector/DetCRD/src/Muon/Muon_Endcap_v01_01.cpp

@@ -269,7 +269,7 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
       superlayer_place.addPhysVolID("Superlayer",i+1);
     }
     dd4hep::PlacedVolume pv;
-    dd4hep::DetElement both_endcap(det_name,2000);
+    dd4hep::DetElement both_endcap(det_name,x_det.id());
     dd4hep::Volume motherVol = theDetector.pickMotherVolume(both_endcap);
     dd4hep::DetElement sdetA = sdet;
     dd4hep::Ref_t(sdetA)->SetName((det_name+"_A").c_str());

Detector/DetCRD/src/Muon/Muon_Endcap_v01_02.cpp

@@ -460,7 +460,7 @@ static dd4hep::Ref_t create_detector(dd4hep::Detector& theDetector,
       superlayer_place.addPhysVolID("Superlayer",i+1);
     }
     dd4hep::PlacedVolume pv;
-    dd4hep::DetElement both_endcap(det_name,2000);
+    dd4hep::DetElement both_endcap(det_name,x_det.id());
     dd4hep::Volume motherVol = theDetector.pickMotherVolume(both_endcap);
     dd4hep::DetElement sdetA = sdet;
     dd4hep::Ref_t(sdetA)->SetName((det_name+"_A").c_str());

Detector/DetIdentifier/include/DetIdentifier/CEPCDetectorData.h

@@ -186,6 +186,41 @@ namespace dd4hep {
 
       return io;
     };
+
+
+    /** Simple data structure with key parameters for
+     *  reconstruction of long crystal bar ECAL
+     *
+     * @author Weizheng Song
+     * @date Nov, 20, 2024
+     * @version $Id: $
+     */
+    struct ECALModuleInfoStruct {
+      int moduleNumber;
+      int staveNumber;
+      int partNumber;
+
+      struct LayerInfo {
+        LayerInfo() :
+          dlayerNumber(-1),
+          slayerNumber(-1),
+          barNumber(-1){
+          }
+        int dlayerNumber;
+        int slayerNumber;
+        int barNumber;
+        };
+      std::vector<LayerInfo> LayerInfos;
+    };
+    typedef StructExtension<ECALModuleInfoStruct> ECALModuleInfoData;
+
+    struct ECALSystemInfoStruct {
+      int systemNumber;
+
+      std::vector<ECALModuleInfoStruct> ModuleInfos;
+    };
+    typedef StructExtension<ECALSystemInfoStruct> ECALSystemInfoData;
+
   }
 }
 

Digitization/DigiCalo/src/EcalDigiAlg.cpp

@@ -109,6 +109,7 @@ StatusCode EcalDigiAlg::initialize()
 		t_SimCont->Branch("step_T1", &m_step_T1);
 		t_SimCont->Branch("step_T2", &m_step_T2);
 		t_SimBar->Branch("totE_Truth", &totE_Truth);
+		t_SimBar->Branch("totE_Truth_MIP", &totE_Truth_MIP);
 		t_SimBar->Branch("totE_Digi", &totE_Digi);
 		t_SimBar->Branch("mean_CT", &mean_CT);
 		t_SimBar->Branch("ECALTemp", &ECALTemp);
@@ -179,6 +180,12 @@ StatusCode EcalDigiAlg::initialize()
 	f_DarkNoise = new TF1("f_DarkNoise", "pow([0]*x, x-1) * exp(-[0]*x) / TMath::Factorial(x)");
 	f_DarkNoise->SetParameter(0, fEcalSiPMCT);
 
+  //ADC non-linearity: a simplified quadratic function, with maximum non-linearity 1%
+  f_ADCNonLin = new TF1("f_ADCNonLin", "pol2", 0, fADC);
+  f_ADCNonLin->SetParameter(0, fADCNonLin/2.);
+  f_ADCNonLin->SetParameter(1, -4*fADCNonLin/fADC);
+  f_ADCNonLin->SetParameter(2, 4*fADCNonLin/fADC/fADC);
+
 	g_CryLYRatio_vs_TID = new TGraph();
 	float TID [8] = {50, 230.043, 731.681, 3448.96, 118953, 1.01164e+7, 1.02341e+8, 2.04907e+8};
 	float CryLYRatio_TID[8] = {0.99, 0.774091, 0.689545, 0.633636, 0.588636, 0.549091, 0.401818, 0.374545};
@@ -220,6 +227,8 @@ StatusCode EcalDigiAlg::execute()
    	if(_nEvt<_Nskip){ _nEvt++; return StatusCode::SUCCESS; }
 
 	Clear();
+  totE_Truth=0;
+  totE_Digi=0;
 
   for(int icol=0; icol<_inputSimHitCollection.size(); icol++){
     try{
@@ -242,7 +251,9 @@ StatusCode EcalDigiAlg::execute()
       
       	if(_Debug>=1) std::cout<<"digi, input sim hit size="<< SimHitCol->size() <<std::endl;
       
+
       totE_Truth=0;
+      totE_Truth_MIP=0;
       totE_Digi=0;
       mean_CT = 0;
       for (int i = 1; i < 10; i++)
@@ -497,10 +508,10 @@ StatusCode EcalDigiAlg::execute()
       //cout<<", SiPM gain "<<sEcalSiPMGainMean<<endl;
       
       		totQ1_Digi = EnergyDigi(ScinGen*totQ1_Att/(totQ1_Att+totQ2_Att), sEcalCryIntLY, sEcalSiPMGainMean, sEcalSiPMDCR, 
-      								f_SiPMResponse, f_SiPMSigmaDet, f_SiPMSigmaRecp, f_SiPMSigmaRecm, f_AsymGauss, f_DarkNoise,
+      								f_SiPMResponse, f_SiPMSigmaDet, f_SiPMSigmaRecp, f_SiPMSigmaRecm, f_AsymGauss, f_DarkNoise, f_ADCNonLin,
       								outLO1, outNDC1, outNDetPE1, outPedestal1, outADC1, outADCGain1)/1000;
       		totQ2_Digi = EnergyDigi(ScinGen*totQ2_Att/(totQ1_Att+totQ2_Att), sEcalCryIntLY, sEcalSiPMGainMean, sEcalSiPMDCR, 
-      								f_SiPMResponse, f_SiPMSigmaDet, f_SiPMSigmaRecp, f_SiPMSigmaRecm, f_AsymGauss, f_DarkNoise,
+      								f_SiPMResponse, f_SiPMSigmaDet, f_SiPMSigmaRecp, f_SiPMSigmaRecm, f_AsymGauss, f_DarkNoise, f_ADCNonLin,
       								outLO2, outNDC2, outNDetPE2, outPedestal2, outADC2, outADCGain2)/1000;
       	}
       	else{
@@ -584,10 +595,11 @@ StatusCode EcalDigiAlg::execute()
       	
         	m_barCol.push_back(hitbar);
       	if(hitbar.getQ1()>=0 && hitbar.getQ2()>=0) totE_Digi += (hitbar.getQ1() + hitbar.getQ2());
-      	//if(totQ1_Att>(fEcalMIPEnergy*fEcalMIP_Thre/1000.) && totQ2_Att>(fEcalMIPEnergy*fEcalMIP_Thre/1000.)){
-      	//  // cout<<"Truth Energy:"<<totQ1_Att<<"   "<<totQ2_Att<<endl;
-      	//  totE_Truth+=(totQ1_Att+totQ2_Att);
-      	//}
+        if (totQ1_Att > (0.001 * fEcalMIPEnergy * fEcalMIP_Thre) && totQ2_Att > (0.001 * fEcalMIPEnergy * fEcalMIP_Thre))
+        {
+//            cout << "Truth Energy:" << totQ1_Att << "   " << totQ2_Att << endl;
+            totE_Truth_MIP += (totQ1_Att + totQ2_Att);
+        }
       
       
       	if(_writeNtuple){
@@ -647,6 +659,7 @@ StatusCode EcalDigiAlg::execute()
 	{
 		std::cout<<"End Loop: Bar Digitalization!"<<std::endl;
 		std::cout<<"Total Truth Energy: "<<totE_Truth<<std::endl;
+        std::cout << "Total Truth Energy with " << (double) fEcalMIP_Thre << " MIP Threshold: " << totE_Truth_MIP << std::endl;
 		std::cout<<"Total Digi Energy: "<<totE_Digi<<std::endl;
 	}
 
@@ -674,12 +687,12 @@ StatusCode EcalDigiAlg::finalize()
 	info() << "Processed " << _nEvt << " events " << endmsg;
   map_readout_decoder.clear();
 	delete m_cellIDConverter, m_geosvc;
-	delete f_SiPMResponse, f_SiPMSigmaDet, f_SiPMSigmaRecp, f_SiPMSigmaRecm, f_AsymGauss, f_DarkNoise, g_SiPMDCR_vs_NIEL, g_CryLYRatio_vs_TID;
+	delete f_SiPMResponse, f_SiPMSigmaDet, f_SiPMSigmaRecp, f_SiPMSigmaRecm, f_AsymGauss, f_DarkNoise, f_ADCNonLin, g_SiPMDCR_vs_NIEL, g_CryLYRatio_vs_TID;
 	return GaudiAlgorithm::finalize();
 }
 
 float EcalDigiAlg::EnergyDigi(float ScinGen, float sEcalCryIntLY, float sEcalSiPMGainMean, float sEcalSiPMDCR, 
-								TF1* f_SiPMResponse, TF1* f_SiPMSigmaDet, TF1* f_SiPMSigmaRecp, TF1* f_SiPMSigmaRecm, TF1* f_AsymGauss, TF1* f_DarkNoise,
+								TF1* f_SiPMResponse, TF1* f_SiPMSigmaDet, TF1* f_SiPMSigmaRecp, TF1* f_SiPMSigmaRecm, TF1* f_AsymGauss, TF1* f_DarkNoise, TF1* f_ADCNonLin, 
 								int& outLO, int& outNDC, int& outNDetPE, float& outPedestal, float& outADC, float& outADCGain)
 {
 	// float sEcalSiPMPDE = rndm.Gaus(fEcalSiPMPDE, fEcalSiPMPDEFlu * fEcalSiPMPDE);
@@ -776,6 +789,7 @@ float EcalDigiAlg::EnergyDigi(float ScinGen, float sEcalCryIntLY, float sEcalSiP
 			float NRec = f_AsymGauss->GetRandom();
 			sADC = NRec * sEcalSiPMGainMean + sPedestal;
 		}
+    sADC = (f_ADCNonLin->Eval(sADC)+1) * sADC;
 		outPedestal = sPedestal;
 		outADCGain = sADC;
 		sPedestal = fPedestal + sEcalSiPMDCR * fEcalTimeInterval * (1 + mean_CT) * sEcalSiPMGainMean;
@@ -820,6 +834,7 @@ float EcalDigiAlg::EnergyDigi(float ScinGen, float sEcalCryIntLY, float sEcalSiP
 			float NRec = f_AsymGauss->GetRandom();
 			sADC = NRec * sEcalSiPMGainMean + sPedestal;
 		}
+    sADC = (f_ADCNonLin->Eval(sADC)+1) * sADC;
 		outPedestal = sPedestal;
 		outADCGain = sADC;
 		sPedestal = fPedestal + sEcalSiPMDCR * fEcalTimeInterval * (1 + mean_CT) * sEcalSiPMGainMean;
@@ -868,6 +883,7 @@ float EcalDigiAlg::EnergyDigi(float ScinGen, float sEcalCryIntLY, float sEcalSiP
 			float NRec = f_AsymGauss->GetRandom();
 			sADC = NRec * sEcalSiPMGainMean + sPedestal;
 		}
+    sADC = (f_ADCNonLin->Eval(sADC)+1) * sADC;
 		outPedestal = sPedestal;
 		outADCGain = sADC;
 		sPedestal = fPedestal + sEcalSiPMDCR * fEcalTimeInterval * (1 + mean_CT) * sEcalSiPMGainMean;
@@ -992,7 +1008,7 @@ edm4hep::MutableSimCalorimeterHit EcalDigiAlg::find(const std::vector<edm4hep::M
 }
 
 void EcalDigiAlg::Clear(){
-  	totE_Truth = -99;
+  totE_Truth = -99;
 	totE_Digi = -99;
 	mean_CT = -99;
 	ECALTemp = -99;

Digitization/DigiCalo/src/EcalDigiAlg.h

@@ -63,7 +63,7 @@ public:
 	edm4hep::MutableSimCalorimeterHit find(const std::vector<edm4hep::MutableSimCalorimeterHit>& m_col, unsigned long long& cellid) const;
   // double Digitization(double edepCry, double fCrosstalkChannel);
   float EnergyDigi(float ScinGen, float sEcalCryMipLY, float sEcalSiPMGainMean, float sEcalSiPMDCR, 
-                    TF1* f_SiPMResponse, TF1* f_SiPMSigmaDet, TF1* f_SiPMSigmaRecp, TF1* f_SiPMSigmaRecm, TF1* f_AsymGauss, TF1* f_DarkNoise,
+                    TF1* f_SiPMResponse, TF1* f_SiPMSigmaDet, TF1* f_SiPMSigmaRecp, TF1* f_SiPMSigmaRecm, TF1* f_AsymGauss, TF1* f_DarkNoise, TF1* f_ADCNonLin,
                     int& outLO, int& outNDC, int& outNDetPE, float& outPedestal, float& outADC, float& outADCGain);
 
 	void Clear();
@@ -81,7 +81,7 @@ protected:
 	TTree* t_SimCont;
 	TTree* t_SimBar;
 
-  double totE_Truth, totE_Digi;
+  double totE_Truth, totE_Truth_MIP, totE_Digi;
   double mean_CT, ECALTemp;
   FloatVec m_step_t;  // yyy: time of each step
 	FloatVec m_step_x, m_step_y, m_step_z, m_step_E, m_step_T1, m_step_T2, m_stepBar_x, m_stepBar_y, m_stepBar_z;
@@ -109,6 +109,7 @@ protected:
   TF1* f_SiPMSigmaRecm = nullptr;
   TF1* f_AsymGauss = nullptr;
   TF1* f_DarkNoise = nullptr;
+  TF1* f_ADCNonLin = nullptr;
   TGraph* g_SiPMDCR_vs_NIEL = nullptr;
   TGraph* g_CryLYRatio_vs_TID = nullptr;
 
@@ -163,6 +164,7 @@ protected:
   mutable Gaudi::Property<int> fSiPMDigiVerbose{this, 	"SiPMDigiVerbose", 1, "SiPM Digitization verbose. 0:w/o response, w/o correction; 1:w/ response, w/o correction; 2:w/ response, w/ simple correction; 3:w/ response, w/ full correction;"};
   mutable Gaudi::Property<float> fGainRatio_12{this, 	"GainRatio_12", 50, "Gain-1 over Gain-2"};
   mutable Gaudi::Property<float> fGainRatio_23{this, 	"GainRatio_23", 60, "Gain-2 over Gain-3"};
+  mutable Gaudi::Property<float> fADCNonLin{this, 	"ADCNonLinearity", 0.01, "ADC non-linearity"};
   mutable Gaudi::Property<float> fEcalMIPEnergy{this, 	"EcalMIPEnergy", 8.9, "MIP energy deposit in 1cm BGO (MeV/MIP)"};
   mutable Gaudi::Property<float> fEcalCryMipLY{this, 	"EcalCryMipLY", 200, "Detected light yield (p.e./MIP)"};
   mutable Gaudi::Property<float> fEcalCryIntLY{this, 	"EcalCryIntLY", 8200, "Intrinsic light yield (ph/MeV)"};