diff --git a/Analysis/MaterialScan/README.md b/Analysis/MaterialScan/README.md
new file mode 100644
index 0000000000000000000000000000000000000000..24cad96167cb7170b57496c92a48e72e9567143b
--- /dev/null
+++ b/Analysis/MaterialScan/README.md
@@ -0,0 +1,84 @@
+# 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 access the interactive interface of materialScan.
+
+```bash
+# Access to the interactive interface of materialScan
+materialScan ./Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01.xml -interactive
+
+# if you want to scan the material for a specific detector, disable/enable the detectors in the xml file
+materialScan ./Detector/DetCRD/compact/TDR_o1_v01/TDR_o1_v01-onlyBeamPipe.xml -interactive
+```
+
+copy the following code to the interactive interface, and press enter.
+*bins* is the number of bins for theta and phi, *world_size* is the size of the scanning area (mm).
+
+```bash
+FILE* outFile = freopen("MaterialScanLogs.txt", "w", stdout);
+int bins = 100;
+double world_size = 10000;
+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);
+ gMaterialScan->print(0,0,0,x,y,z);
+ }
+}
+```
+
+The "FILE* outFile" truncates the "stdout" output stream to the file "MaterialScanLogs.txt", 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.
+
+The file "MaterialScanLogs.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
diff --git a/Analysis/MaterialScan/options/extract_multi_x0.py b/Analysis/MaterialScan/options/extract_multi_x0.py
new file mode 100644
index 0000000000000000000000000000000000000000..b8e685c15f2f9feb2ee0904b76e19a5e2a33c2d5
--- /dev/null
+++ b/Analysis/MaterialScan/options/extract_multi_x0.py
@@ -0,0 +1,116 @@
+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
diff --git a/Analysis/MaterialScan/options/extract_x0.py b/Analysis/MaterialScan/options/extract_x0.py
new file mode 100644
index 0000000000000000000000000000000000000000..d53f9350cf8999e4dc989b3608080e42ff9650c4
--- /dev/null
+++ b/Analysis/MaterialScan/options/extract_x0.py
@@ -0,0 +1,94 @@
+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