From 12cbfa02648d2ded184e6f878b6572aceb71a79c Mon Sep 17 00:00:00 2001 From: kkaan Date: Fri, 7 Jun 2024 12:40:46 +1000 Subject: [PATCH 1/8] initial guess to prevent overflow errors in curve_fit --- requirements.txt | 9 ++- .../correction_coefficients.py | 60 ++++++++++++------- 2 files changed, 45 insertions(+), 24 deletions(-) diff --git a/requirements.txt b/requirements.txt index 17db2d2..a6f6dec 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,5 +1,8 @@ sphinx sphinx-rtd-theme -numpy -pandas -seaborn +numpy~=1.26.4 +pandas~=2.2.1 +seaborn~=0.13.2 + +matplotlib~=3.8.4 +scipy~=1.13.1 \ No newline at end of file diff --git a/src/correction_coefficients/correction_coefficients.py b/src/correction_coefficients/correction_coefficients.py index bce771e..8cb729c 100644 --- a/src/correction_coefficients/correction_coefficients.py +++ b/src/correction_coefficients/correction_coefficients.py @@ -4,7 +4,7 @@ The module includes the following functions: - `get_pr_data_from_acm_files`: Extracts the average counts/50ms from ACM files in a specified directory. -- `exponential_fit`: Defines an exponential fitting function for curve fitting. +- `saturation_fit`: Defines an exponential fitting function for curve fitting. - `get_correction_coefficients`: Calculates the correction coefficients based on the average counts/50ms and relative signal values. - `plot_correction_curve`: Plots the correction curve with the data points and the fitted exponential function. - `plot_counts_per_50ms`: Plots the counts per 50ms at different nominal dose rates. @@ -17,6 +17,7 @@ import os import sys +import pandas as pd import numpy as np import matplotlib.pyplot as plt import seaborn as sns @@ -42,36 +43,50 @@ def get_pr_data_from_acm_files(pr_measurement_folder): List of tuples containing the file name and the average counts/50ms for each nominal dose rate. """ diode_numbers = np.array([758, 759, 760, 761, 692, 693, 694, 695, 626, 627, 628, 629]).astype(str) + + # initialize a dictionary to store the average counts/50ms for each file over all + # diodes in the middle 200 frames. + counts_per_50ms_dict = {} + + # initialize a list to store the average counts per 50ms for each file counts_per_50ms_list = [] for acm_file in os.listdir(pr_measurement_folder): acm_file_path = os.path.join(pr_measurement_folder, acm_file) frame_data_df, diode_data_df, bkrnd_and_calibration_df = io_snc.parse_acm_file(acm_file_path) - # Ensure the frame data is sorted by time or frame number if necessary - diode_data_df.sort_index(inplace=True) - - # Calculate the count rate by differencing the accumulated counts - diode_data_diff = diode_data_df.diff().iloc[1:] - # Grab the count rate from the middle 100 frames of the measurement - num_frames = len(diode_data_diff) + num_frames = len(diode_data_df) start_frame = num_frames // 2 - 50 end_frame = start_frame + 100 - middle_frames = diode_data_diff.iloc[start_frame:end_frame] + middle_frames = diode_data_df.iloc[start_frame:end_frame] + + # Calculate the count rate by differencing the accumulated counts + diode_data_diff_df = middle_frames.diff().iloc[1:] # Extract the count rate for the specified diodes in the middle frames - count_rate_df = middle_frames.loc[:, diode_numbers] + count_rate_df = diode_data_diff_df.loc[:, diode_numbers] + + # Calculate the average count rate for each frame across the selected diodes + avg_count_rate_per_frame = count_rate_df.mean(axis=1) + avg_count_rate_per_frame = avg_count_rate_per_frame.values + counts_per_50ms_dict[acm_file] = avg_count_rate_per_frame # Calculate the average count rate for the selected diodes avg_count_rate = count_rate_df.mean().mean() counts_per_50ms_list.append((acm_file, avg_count_rate)) + + # Convert the dictionary to a DataFrame + counts_per_50ms_df = pd.DataFrame.from_dict(counts_per_50ms_dict, orient='index') + # Save the DataFrame to a CSV file + counts_per_50ms_df.to_csv('counts_per_50ms.csv') + counts_per_50ms_list = sorted(counts_per_50ms_list, key=lambda x: x[1]) return counts_per_50ms_list -def exponential_fit(x, a, b, c): +def saturation_func(x, a, b, c): """ Exponential fitting function. @@ -91,7 +106,7 @@ def exponential_fit(x, a, b, c): array_like The calculated dependent variable values. """ - return a * np.exp(b * x) + c + return c - a * np.exp(-b * x) def get_correction_coefficients(counts_per_50ms, relative_signal): @@ -111,10 +126,12 @@ def get_correction_coefficients(counts_per_50ms, relative_signal): The fitted coefficients (a, b, c) for the exponential correction function. """ counts_per_50ms = np.array([data[1] for data in counts_per_50ms]) + initial_guess = [0.035, 5.21e-5, 1.0] + params, covariance = curve_fit(saturation_func, counts_per_50ms, + relative_signal, p0=initial_guess) + a_fit, b_fit, c_fit = params - # Perform the curve fitting - popt, _ = curve_fit(exponential_fit, counts_per_50ms, relative_signal) - return tuple(popt) + return a_fit, b_fit, c_fit, covariance def plot_correction_curve(counts_per_50ms, relative_signal, a_fit, b_fit, c_fit): @@ -138,7 +155,7 @@ def plot_correction_curve(counts_per_50ms, relative_signal, a_fit, b_fit, c_fit) # Generate fitted curve data points x_fit = np.linspace(counts_per_50ms.min(), counts_per_50ms.max(), 500) - y_fit = exponential_fit(x_fit, a_fit, b_fit, c_fit) + y_fit = saturation_func(x_fit, a_fit, b_fit, c_fit) # Set Seaborn style sns.set(style="whitegrid") @@ -162,8 +179,8 @@ def plot_correction_curve(counts_per_50ms, relative_signal, a_fit, b_fit, c_fit) plt.title('Pulse rate (MU/min)') plt.legend() plt.grid(True) - plt.show() + plt.savefig('counts_per_50ms.png') def plot_counts_per_50ms(counts_per_50ms): @@ -179,7 +196,7 @@ def plot_counts_per_50ms(counts_per_50ms): # Set Seaborn style sns.set(style="whitegrid") - + plt.ioff() # Turn off interactive mode # Create a plot with Seaborn plt.figure(figsize=(10, 6)) sns.scatterplot(x=range(len(counts_per_50ms)), y=counts_per_50ms, label='Counts/50ms') @@ -190,8 +207,9 @@ def plot_counts_per_50ms(counts_per_50ms): plt.title('Counts per 50ms at different nominal dose rates') plt.legend() plt.grid(True) - plt.show() + plt.savefig('counts_per_50ms.png') + def main(): @@ -200,9 +218,9 @@ def main(): """ pr_measurement_folder = r"P:\02_QA Equipment\02_ArcCheck\05_Commissoning\03_NROAC\Dose Rate Dependence Fix\NRO PR Coefficient" counts_per_50ms = get_pr_data_from_acm_files(pr_measurement_folder) - plot_counts_per_50ms(counts_per_50ms) + # plot_counts_per_50ms(counts_per_50ms) relative_signal = np.array([0.964, 0.971, 0.980, 0.988, 0.994, 0.995, 1.000, 1.000]) - a_fit, b_fit, c_fit = get_correction_coefficients(counts_per_50ms, relative_signal) + a_fit, b_fit, c_fit, covariance = get_correction_coefficients(counts_per_50ms, relative_signal) plot_correction_curve(counts_per_50ms, relative_signal, a_fit, b_fit, c_fit) From 23a4275c9de20cbe8debe77d8e9537398260b240 Mon Sep 17 00:00:00 2001 From: kkaan Date: Fri, 7 Jun 2024 12:41:35 +1000 Subject: [PATCH 2/8] adding some placeholder files for utility ideas --- src/Utilities/_init_.py | 0 src/Utilities/fix_no_dose_error_composite.py | 22 ++++++++++++++++++++ 2 files changed, 22 insertions(+) create mode 100644 src/Utilities/_init_.py create mode 100644 src/Utilities/fix_no_dose_error_composite.py diff --git a/src/Utilities/_init_.py b/src/Utilities/_init_.py new file mode 100644 index 0000000..e69de29 diff --git a/src/Utilities/fix_no_dose_error_composite.py b/src/Utilities/fix_no_dose_error_composite.py new file mode 100644 index 0000000..e22039d --- /dev/null +++ b/src/Utilities/fix_no_dose_error_composite.py @@ -0,0 +1,22 @@ +import sys +import os + +# Add the src directory to the sys.path +sys.path.append(os.path.abspath(os.path.join(os.path.dirname(__file__), '..'))) +from src import io_snc + +def main(): + """ + main function to fix measured files + + Returns + ------- + + """ + # ask user for file path of measured file + file_path = input("Please enter the file path: ") + return file_path + + +if __name__ == "__main__": + main() \ No newline at end of file From 0859bf2c4c0d262e2675fe8c4466f24fb9202e21 Mon Sep 17 00:00:00 2001 From: kkaan Date: Wed, 12 Jun 2024 14:07:10 +1000 Subject: [PATCH 3/8] plottings --- src/plots_for_presentation/difference_plot.py | 162 ++++++++++++++++++ 1 file changed, 162 insertions(+) create mode 100644 src/plots_for_presentation/difference_plot.py diff --git a/src/plots_for_presentation/difference_plot.py b/src/plots_for_presentation/difference_plot.py new file mode 100644 index 0000000..425c446 --- /dev/null +++ b/src/plots_for_presentation/difference_plot.py @@ -0,0 +1,162 @@ +import pandas as pd +import seaborn as sns +import matplotlib.pyplot as plt +# Reduce the font size + +# Load data +file_path = r'P:\02_QA Equipment\02_ArcCheck\05_Commissoning\03_NROAC\Dose Rate Dependence Fix\Test on script\ACResultsPowerQuerySummary.txt' +df = pd.read_csv(file_path, delimiter='\t') + +# Ensure numeric columns are converted from strings to numbers +for col in ['3%3mm', '2%3mm', '1%3mm']: + df[col] = pd.to_numeric(df[col], errors='coerce') + +# Split data into Original and PR correction +original_df = df[df['Correction'] == 'Original'].set_index('Plan Name') +pr_df = df[df['Correction'] == 'PR'].set_index('Plan Name') + + +# Check if data is numeric +if (original_df[['3%3mm', '2%3mm', '1%3mm']].dtypes == 'float64').all() and (pr_df[['3%3mm', '2%3mm', '1%3mm']].dtypes == 'float64').all(): + # Calculate differences + difference_df = pr_df[['3%3mm', '2%3mm', '1%3mm']].subtract(original_df[['3%3mm', '2%3mm', '1%3mm']]) +else: + print("Data is not numeric. Cannot perform subtraction.") + + +# Reset index to use for plotting +difference_df.reset_index(inplace=True) + +# Melt the dataframe for seaborn plotting +melted_df = difference_df.melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='Difference') + +# Define a function to map the improvement/deterioration to colors +def color_mapper(val): + if val >= 0: + return 'green' + else: + return 'red' + +# Melt the original and PR dataframes +original_melted = original_df.reset_index().melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='Original') +pr_melted = pr_df.reset_index().melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='PR') + +import matplotlib as mpl + +import matplotlib as mpl + +# Create the plot +fig, axs = plt.subplots(3, 1, figsize=(10, 24)) +sns.set(style="whitegrid") + +metrics = ['3%3mm', '2%3mm', '1%3mm'] + +for i, metric in enumerate(metrics): + ax1 = axs[i] + + # Plot the differences + subset_df = melted_df[melted_df['Metric'] == metric] + colors = subset_df['Difference'].apply(color_mapper) + sns.scatterplot(x=subset_df['Difference'], y=subset_df['Plan Name'], hue=colors, palette=['red', 'green'], legend=False, s=100, marker='o', ax=ax1) + + # Set the limits for the first x-axis and set label color to dark grey + ax1.set_xlim([-10, 10]) + ax1.set_xlabel('Differences', color='darkgrey') + + # Reduce the font size + mpl.rcParams.update({'font.size': 10}) + + # Create a second x-axis + ax2 = ax1.twiny() + + # Plot the original values + subset_df = original_melted[original_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['Original'], y=subset_df['Plan Name'], color='darkgrey', legend=False, s=100, marker='o', alpha=0.5, ax=ax2) + + # Plot the PR values + subset_df = pr_melted[pr_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['PR'], y=subset_df['Plan Name'], color='darkblue', legend=False, s=100, marker='o', alpha=0.5, ax=ax2) + + # Set the limits for the second x-axis and set label color to dark grey + ax2.set_xlim([50, 110]) + ax2.set_xlabel('Pass Rates', color='darkgrey') + + # Add labels and title + ax1.set_ylabel('Plan Name') + ax1.set_title(metric, loc='left') # Set title to just the metric and align left + ax1.axvline(0, color='gray', linestyle='--') + +# Adjust the space between the plots +plt.subplots_adjust(hspace=0.5) + +plt.tight_layout() +plt.ioff() +plt.savefig('pr_correction_plot.png', dpi=300) +plt.show() + + +# Load data for DPP correction +dpp_df = df[df['Correction'] == 'DPP'].set_index('Plan Name') + +# Check if data is numeric +if (original_df[['3%3mm', '2%3mm', '1%3mm']].dtypes == 'float64').all() and (dpp_df[['3%3mm', '2%3mm', '1%3mm']].dtypes == 'float64').all(): + # Calculate differences + difference_dpp_df = dpp_df[['3%3mm', '2%3mm', '1%3mm']].subtract(original_df[['3%3mm', '2%3mm', '1%3mm']]) +else: + print("Data is not numeric. Cannot perform subtraction.") + +# Reset index to use for plotting +difference_dpp_df.reset_index(inplace=True) + +# Melt the dataframe for seaborn plotting +melted_dpp_df = difference_dpp_df.melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='Difference') + +# Melt the DPP dataframe +dpp_melted = dpp_df.reset_index().melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='DPP') + +# Create the plot for DPP +fig, axs = plt.subplots(3, 1, figsize=(10, 24)) +sns.set(style="whitegrid") + +for i, metric in enumerate(metrics): + ax1 = axs[i] + + # Plot the differences + subset_df = melted_dpp_df[melted_dpp_df['Metric'] == metric] + colors = subset_df['Difference'].apply(color_mapper) + sns.scatterplot(x=subset_df['Difference'], y=subset_df['Plan Name'], hue=colors, palette=['green', 'red'], legend=False, s=100, marker='o', ax=ax1) + + # Set the limits for the first x-axis and set label color to dark grey + ax1.set_xlim([-10, 10]) + ax1.set_xlabel('Differences', color='darkgrey') + + # Reduce the font size + mpl.rcParams.update({'font.size': 10}) + + # Create a second x-axis + ax2 = ax1.twiny() + + # Plot the original values + subset_df = original_melted[original_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['Original'], y=subset_df['Plan Name'], color='darkgrey', legend=False, s=100, marker='o', alpha=0.5, ax=ax2) + + # Plot the DPP values + subset_df = dpp_melted[dpp_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['DPP'], y=subset_df['Plan Name'], color='darkblue', legend=False, s=100, marker='o', alpha=0.5, ax=ax2) + + # Set the limits for the second x-axis and set label color to dark grey + ax2.set_xlim([50, 110]) + ax2.set_xlabel('Pass Rates', color='darkgrey') + + # Add labels and title + ax1.set_ylabel('Plan Name') + ax1.set_title(metric, loc='left') # Set title to just the metric and align left + ax1.axvline(0, color='gray', linestyle='--') + +# Adjust the space between the plots +plt.subplots_adjust(hspace=0.5) + +plt.tight_layout() +plt.ioff() +plt.savefig('dpp_correction_plot.png', dpi=300) +plt.show() \ No newline at end of file From 119a26e426d843bf83134e8e7e0799cb74d0869f Mon Sep 17 00:00:00 2001 From: kkaan Date: Wed, 12 Jun 2024 14:38:19 +1000 Subject: [PATCH 4/8] plottings for pdfing --- src/plots_for_presentation/difference_plot.py | 87 ++++++++++--------- 1 file changed, 46 insertions(+), 41 deletions(-) diff --git a/src/plots_for_presentation/difference_plot.py b/src/plots_for_presentation/difference_plot.py index 425c446..74cb647 100644 --- a/src/plots_for_presentation/difference_plot.py +++ b/src/plots_for_presentation/difference_plot.py @@ -1,7 +1,7 @@ import pandas as pd import seaborn as sns import matplotlib.pyplot as plt -# Reduce the font size +import matplotlib as mpl # Load data file_path = r'P:\02_QA Equipment\02_ArcCheck\05_Commissoning\03_NROAC\Dose Rate Dependence Fix\Test on script\ACResultsPowerQuerySummary.txt' @@ -52,12 +52,14 @@ def color_mapper(val): metrics = ['3%3mm', '2%3mm', '1%3mm'] for i, metric in enumerate(metrics): - ax1 = axs[i] + # Create the plot for each metric + fig, ax1 = plt.subplots(figsize=(10, 8)) + sns.set(style="whitegrid") # Plot the differences subset_df = melted_df[melted_df['Metric'] == metric] colors = subset_df['Difference'].apply(color_mapper) - sns.scatterplot(x=subset_df['Difference'], y=subset_df['Plan Name'], hue=colors, palette=['red', 'green'], legend=False, s=100, marker='o', ax=ax1) + sns.scatterplot(x=subset_df['Difference'], y=subset_df['Plan Name'], hue=colors, palette=['red', 'green'], legend=False, s=100, marker='^', ax=ax1) # Set the limits for the first x-axis and set label color to dark grey ax1.set_xlim([-10, 10]) @@ -83,16 +85,14 @@ def color_mapper(val): # Add labels and title ax1.set_ylabel('Plan Name') - ax1.set_title(metric, loc='left') # Set title to just the metric and align left + ax1.set_title(f'PR {metric}metric', loc='left') # Set title to just the metric and align left ax1.axvline(0, color='gray', linestyle='--') -# Adjust the space between the plots -plt.subplots_adjust(hspace=0.5) + plt.tight_layout() + plt.ioff() + plt.savefig(f'PR_{metric}_correction_plot.png', dpi=300) + -plt.tight_layout() -plt.ioff() -plt.savefig('pr_correction_plot.png', dpi=300) -plt.show() # Load data for DPP correction @@ -119,44 +119,49 @@ def color_mapper(val): sns.set(style="whitegrid") for i, metric in enumerate(metrics): - ax1 = axs[i] + metrics = ['3%3mm', '2%3mm', '1%3mm'] - # Plot the differences - subset_df = melted_dpp_df[melted_dpp_df['Metric'] == metric] - colors = subset_df['Difference'].apply(color_mapper) - sns.scatterplot(x=subset_df['Difference'], y=subset_df['Plan Name'], hue=colors, palette=['green', 'red'], legend=False, s=100, marker='o', ax=ax1) + for i, metric in enumerate(metrics): + # Create the plot for each metric + fig, ax1 = plt.subplots(figsize=(10, 8)) + sns.set(style="whitegrid") - # Set the limits for the first x-axis and set label color to dark grey - ax1.set_xlim([-10, 10]) - ax1.set_xlabel('Differences', color='darkgrey') + # Plot the differences + subset_df = melted_dpp_df[melted_df['Metric'] == metric] + colors = subset_df['Difference'].apply(color_mapper) + sns.scatterplot(x=subset_df['Difference'], y=subset_df['Plan Name'], hue=colors, palette=['red', 'green'], + legend=False, s=100, marker='^', ax=ax1) - # Reduce the font size - mpl.rcParams.update({'font.size': 10}) + # Set the limits for the first x-axis and set label color to dark grey + ax1.set_xlim([-10, 10]) + ax1.set_xlabel('Differences', color='darkgrey') - # Create a second x-axis - ax2 = ax1.twiny() + # Reduce the font size + mpl.rcParams.update({'font.size': 10}) - # Plot the original values - subset_df = original_melted[original_melted['Metric'] == metric] - sns.scatterplot(x=subset_df['Original'], y=subset_df['Plan Name'], color='darkgrey', legend=False, s=100, marker='o', alpha=0.5, ax=ax2) + # Create a second x-axis + ax2 = ax1.twiny() - # Plot the DPP values - subset_df = dpp_melted[dpp_melted['Metric'] == metric] - sns.scatterplot(x=subset_df['DPP'], y=subset_df['Plan Name'], color='darkblue', legend=False, s=100, marker='o', alpha=0.5, ax=ax2) + # Plot the original values + subset_df = original_melted[original_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['Original'], y=subset_df['Plan Name'], color='darkgrey', legend=False, s=100, + marker='o', alpha=0.5, ax=ax2) - # Set the limits for the second x-axis and set label color to dark grey - ax2.set_xlim([50, 110]) - ax2.set_xlabel('Pass Rates', color='darkgrey') + # Plot the PR values + subset_df = dpp_melted[dpp_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['DPP'], y=subset_df['Plan Name'], color='darkblue', legend=False, s=100, marker='o', + alpha=0.5, ax=ax2) - # Add labels and title - ax1.set_ylabel('Plan Name') - ax1.set_title(metric, loc='left') # Set title to just the metric and align left - ax1.axvline(0, color='gray', linestyle='--') + # Set the limits for the second x-axis and set label color to dark grey + ax2.set_xlim([50, 110]) + ax2.set_xlabel('Pass Rates', color='darkgrey') + + # Add labels and title + ax1.set_ylabel('Plan Name') + ax1.set_title(f'DPP {metric}', loc='left') # Set title to just the metric and align left + ax1.axvline(0, color='gray', linestyle='--') -# Adjust the space between the plots -plt.subplots_adjust(hspace=0.5) + plt.tight_layout() + plt.ioff() + plt.savefig(f'DPP_{metric}_correction_plot.png', dpi=300) -plt.tight_layout() -plt.ioff() -plt.savefig('dpp_correction_plot.png', dpi=300) -plt.show() \ No newline at end of file From 1ee1199ac26ee4cba58da9e24fe836d7ab1a2cf9 Mon Sep 17 00:00:00 2001 From: kkaan Date: Wed, 12 Jun 2024 14:38:38 +1000 Subject: [PATCH 5/8] correction coefficients extracted from measurement --- .../correction_coefficients.py | 63 +++++++++++++++---- src/corrections.py | 4 +- src/main.py | 1 + src/plots_for_presentation/summaryboxplots.py | 12 ++-- 4 files changed, 64 insertions(+), 16 deletions(-) diff --git a/src/correction_coefficients/correction_coefficients.py b/src/correction_coefficients/correction_coefficients.py index 8cb729c..635ba87 100644 --- a/src/correction_coefficients/correction_coefficients.py +++ b/src/correction_coefficients/correction_coefficients.py @@ -17,11 +17,14 @@ import os import sys -import pandas as pd -import numpy as np import matplotlib.pyplot as plt +import numpy as np +import pandas as pd import seaborn as sns from scipy.optimize import curve_fit +import datetime + + # Add the src directory to the sys.path sys.path.append(os.path.abspath(os.path.join(os.path.dirname(__file__), '..'))) @@ -76,7 +79,6 @@ def get_pr_data_from_acm_files(pr_measurement_folder): avg_count_rate = count_rate_df.mean().mean() counts_per_50ms_list.append((acm_file, avg_count_rate)) - # Convert the dictionary to a DataFrame counts_per_50ms_df = pd.DataFrame.from_dict(counts_per_50ms_dict, orient='index') # Save the DataFrame to a CSV file @@ -127,12 +129,49 @@ def get_correction_coefficients(counts_per_50ms, relative_signal): """ counts_per_50ms = np.array([data[1] for data in counts_per_50ms]) initial_guess = [0.035, 5.21e-5, 1.0] - params, covariance = curve_fit(saturation_func, counts_per_50ms, - relative_signal, p0=initial_guess) + # noinspection PyTupleAssignmentBalance + params, covariance = curve_fit(saturation_func, counts_per_50ms, + relative_signal, p0=initial_guess) a_fit, b_fit, c_fit = params return a_fit, b_fit, c_fit, covariance +def get_coefficient_file_path(date=None): + """ + Get the path to the file where the correction coefficients will be written. + + Returns + ------- + str + The file path for the correction coefficients. + """ + # get user input on type of whether it DPP or PR + coefficient_file_directory = input("Enter the directory to the file where the correction coefficients will be written: ") + + # Get the DRD type, ArcCheck SN from the user + DRD_type = input("Enter the DRD type: ") + ArcCheck_SN = input("Enter the ArcCheck SN: ") + + # Get today's date + today = datetime.date.today() + date = today.strftime("%Y-%m-%d") + + # construct coefficient file name from DRD type, ArcCheck SN and Date + coefficient_file_name = f"{DRD_type}_{ArcCheck_SN}_{date}.txt" + + # construct the full path to the coefficient file + coefficient_file_path = os.path.join(coefficient_file_directory, coefficient_file_name) + + return coefficient_file_path + +def write_correction_coefficients_to_file(coefficient_file_path, a_fit, b_fit, c_fit, covariance): + with open(coefficient_file_path, 'w') as f: + # Write the variable names and their values to the file + f.write(f"a_fit = {a_fit}\n") + f.write(f"b_fit = {b_fit}\n") + f.write(f"c_fit = {c_fit}\n") + f.write(f"covariance = {covariance}\n") + def plot_correction_curve(counts_per_50ms, relative_signal, a_fit, b_fit, c_fit): """ @@ -165,7 +204,7 @@ def plot_correction_curve(counts_per_50ms, relative_signal, a_fit, b_fit, c_fit) sns.lineplot(x=x_fit, y=y_fit, label=f'Exponential fit, a={a_fit:.4f}, b={b_fit:.2e}, c={c_fit:.3f}') sns.scatterplot(x=counts_per_50ms, y=relative_signal, label='Data points') - plt.ioff() # Turn off interactive mode + #plt.ioff() # Turn off interactive mode # Plot the data points and the fitted curve plt.figure(figsize=(8, 6)) @@ -196,7 +235,7 @@ def plot_counts_per_50ms(counts_per_50ms): # Set Seaborn style sns.set(style="whitegrid") - plt.ioff() # Turn off interactive mode + #plt.ioff() # Turn off interactive mode # Create a plot with Seaborn plt.figure(figsize=(10, 6)) sns.scatterplot(x=range(len(counts_per_50ms)), y=counts_per_50ms, label='Counts/50ms') @@ -208,20 +247,22 @@ def plot_counts_per_50ms(counts_per_50ms): plt.legend() plt.grid(True) plt.show() - plt.savefig('counts_per_50ms.png') - + plt.savefig('counts_per_50ms_nominal.png') def main(): """ Main function to execute the correction coefficient calculations and plotting. """ + pr_measurement_folder = r"P:\02_QA Equipment\02_ArcCheck\05_Commissoning\03_NROAC\Dose Rate Dependence Fix\NRO PR Coefficient" counts_per_50ms = get_pr_data_from_acm_files(pr_measurement_folder) # plot_counts_per_50ms(counts_per_50ms) relative_signal = np.array([0.964, 0.971, 0.980, 0.988, 0.994, 0.995, 1.000, 1.000]) - a_fit, b_fit, c_fit, covariance = get_correction_coefficients(counts_per_50ms, relative_signal) - plot_correction_curve(counts_per_50ms, relative_signal, a_fit, b_fit, c_fit) + a_fit, b_fit, c_fit, covariance = get_correction_coefficients(counts_per_50ms, relative_signal) + coefficient_file_path = get_coefficient_file_path() + write_correction_coefficients_to_file(coefficient_file_path, a_fit, b_fit, c_fit, covariance) + # plot_correction_curve(counts_per_50ms, relative_signal, a_fit, b_fit, c_fit) if __name__ == "__main__": diff --git a/src/corrections.py b/src/corrections.py index 1ff9adf..73a6703 100644 --- a/src/corrections.py +++ b/src/corrections.py @@ -73,7 +73,9 @@ def get_intrinsic_corrections(array_data): def apply_jager_corrections(counts_accumulated_df, bkrnd_and_calibration_df, intrinsic_corrections=None): """Apply Jager pulse rate and dose per pulse corrections.""" - a_pr, b_pr, c_pr = 0.035, 5.21 * 10 ** -5, 1 + + # a_pr, b_pr, c_pr = 0.035, 5.21 * 10 ** -5, 1 # Jager coefficients from paper + a_pr, b_pr, c_pr = 0.03921672747199777, 6.499040838838414e-055, 0.999142799734291 # Coefficients from NRO Measurements jager_pr_coefficients = np.array([a_pr, b_pr, c_pr]) a_dpp, b_dpp, c_dpp = 0.0978, 3.33 * 10 ** -5, 1.011 diff --git a/src/main.py b/src/main.py index 727a98e..1dd16fa 100644 --- a/src/main.py +++ b/src/main.py @@ -199,3 +199,4 @@ def main(): if __name__ == "__main__": main() + diff --git a/src/plots_for_presentation/summaryboxplots.py b/src/plots_for_presentation/summaryboxplots.py index ece3faa..ecb73d0 100644 --- a/src/plots_for_presentation/summaryboxplots.py +++ b/src/plots_for_presentation/summaryboxplots.py @@ -7,23 +7,27 @@ data = pd.read_csv(file_path, sep='\t') # Melt the data -data_melted = data.melt(id_vars=['Plan Name', '(mm)', 'Custom'], +data_melted = data.melt(id_vars=['Plan Name', '(mm)', 'Correction'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Percentage', value_name='Value') +# Filter the data to only include 'Original' and 'DPP' +data_melted_filtered = data_melted[data_melted['Correction'].isin(['Original', 'DPP'])] + # Set the theme sns.set_theme(style="white") # Set up the matplotlib figure with smaller size plt.figure(figsize=(20, 5)) +plt.figure(figsize=(20, 5)) # Creating box plots with a custom color palette -palette = {'Original': '#1f77b4', 'DPP': '#2ca02c', 'PR': '#17becf', 'PRI': '#ff7f0e'} -y_limits = (75, data_melted['Value'].max() + 5) # Set minimum value to 75 and adjust the max limit for clarity +palette = {'Original': '#1f77b4', 'DPP': '#2ca02c'} +y_limits = (75, data_melted_filtered['Value'].max() + 5) # Set minimum value to 75 and adjust the max limit for clarity for i, percentage in enumerate(['3%3mm', '2%3mm', '1%3mm'], start=1): plt.subplot(1, 3, i) - sns.boxplot(x='Custom', y='Value', data=data_melted[data_melted['Percentage'] == percentage], palette=palette) + sns.boxplot(x='Correction', y='Value', data=data_melted_filtered[data_melted_filtered['Percentage'] == percentage], palette=palette) plt.title(f'{percentage}', fontsize=12) # Increase title font size plt.ylim(y_limits) plt.xticks(rotation=45, fontsize=10) # Increase x-tick labels font size From f26d8ca579dda735ca1a7a58860d85a39cca1007 Mon Sep 17 00:00:00 2001 From: kkaan Date: Mon, 17 Jun 2024 13:29:20 +1000 Subject: [PATCH 6/8] Plotting for presentation --- src/plots_for_presentation/difference_plot.py | 4 +- .../difference_plot_dpp.py | 46 +++++++++++++++++++ .../difference_plot_pr.py | 46 +++++++++++++++++++ src/plots_for_presentation/legend.py | 19 ++++++++ 4 files changed, 113 insertions(+), 2 deletions(-) create mode 100644 src/plots_for_presentation/difference_plot_dpp.py create mode 100644 src/plots_for_presentation/difference_plot_pr.py create mode 100644 src/plots_for_presentation/legend.py diff --git a/src/plots_for_presentation/difference_plot.py b/src/plots_for_presentation/difference_plot.py index 74cb647..d545344 100644 --- a/src/plots_for_presentation/difference_plot.py +++ b/src/plots_for_presentation/difference_plot.py @@ -85,7 +85,7 @@ def color_mapper(val): # Add labels and title ax1.set_ylabel('Plan Name') - ax1.set_title(f'PR {metric}metric', loc='left') # Set title to just the metric and align left + ax1.set_title(f'PR {metric}', loc='left') # Set title to just the metric and align left ax1.axvline(0, color='gray', linestyle='--') plt.tight_layout() @@ -129,7 +129,7 @@ def color_mapper(val): # Plot the differences subset_df = melted_dpp_df[melted_df['Metric'] == metric] colors = subset_df['Difference'].apply(color_mapper) - sns.scatterplot(x=subset_df['Difference'], y=subset_df['Plan Name'], hue=colors, palette=['red', 'green'], + sns.scatterplot(x=subset_df['Difference'], y=subset_df['Plan Name'], hue=colors, palette=['green', 'red'], legend=False, s=100, marker='^', ax=ax1) # Set the limits for the first x-axis and set label color to dark grey diff --git a/src/plots_for_presentation/difference_plot_dpp.py b/src/plots_for_presentation/difference_plot_dpp.py new file mode 100644 index 0000000..e2bb8b2 --- /dev/null +++ b/src/plots_for_presentation/difference_plot_dpp.py @@ -0,0 +1,46 @@ +import pandas as pd +import seaborn as sns +import matplotlib.pyplot as plt + +# Load data +file_path = r'P:\02_QA Equipment\02_ArcCheck\05_Commissoning\03_NROAC\Dose Rate Dependence Fix\Test on script\ACResultsPowerQuerySummary.txt' +df = pd.read_csv(file_path, delimiter='\t') + +# Ensure numeric columns are converted from strings to numbers +for col in ['3%3mm', '2%3mm', '1%3mm']: + df[col] = pd.to_numeric(df[col], errors='coerce') + +# Split data into Original and DPP correction +original_df = df[df['Correction'] == 'Original'].set_index('Plan Name') +dpp_df = df[df['Correction'] == 'DPP'].set_index('Plan Name') + +# Melt the original and DPP dataframes +original_melted = original_df.reset_index().melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='Original') +dpp_melted = dpp_df.reset_index().melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='DPP') + +metrics = ['3%3mm', '2%3mm', '1%3mm'] + +for i, metric in enumerate(metrics): + # Create the plot for each metric + fig, ax = plt.subplots(figsize=(10, 8)) + sns.set(style="whitegrid") + + # Plot the original values + subset_df = original_melted[original_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['Original'], y=subset_df['Plan Name'], color='darkgrey', legend=False, s=100, marker='o', alpha=0.5, ax=ax) + + # Plot the DPP values + subset_df = dpp_melted[dpp_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['DPP'], y=subset_df['Plan Name'], color='darkblue', legend=False, s=100, marker='o', alpha=0.5, ax=ax) + + # Set the limits for the x-axis and set label color to dark grey + ax.set_xlim([50, 110]) + ax.set_xlabel('Pass Rates', color='darkgrey') + + # Add labels and title + ax.set_ylabel('Plan Name') + ax.set_title(f'DPP {metric}', loc='left') # Set title to just the metric and align left + + plt.tight_layout() + plt.savefig(f'DPP_{metric}_pass_rate_plot.png', dpi=300) + plt.show() \ No newline at end of file diff --git a/src/plots_for_presentation/difference_plot_pr.py b/src/plots_for_presentation/difference_plot_pr.py new file mode 100644 index 0000000..0ab42e7 --- /dev/null +++ b/src/plots_for_presentation/difference_plot_pr.py @@ -0,0 +1,46 @@ +import pandas as pd +import seaborn as sns +import matplotlib.pyplot as plt + +# Load data +file_path = r'P:\02_QA Equipment\02_ArcCheck\05_Commissoning\03_NROAC\Dose Rate Dependence Fix\Test on script\ACResultsPowerQuerySummary.txt' +df = pd.read_csv(file_path, delimiter='\t') + +# Ensure numeric columns are converted from strings to numbers +for col in ['3%3mm', '2%3mm', '1%3mm']: + df[col] = pd.to_numeric(df[col], errors='coerce') + +# Split data into Original and PR correction +original_df = df[df['Correction'] == 'Original'].set_index('Plan Name') +pr_df = df[df['Correction'] == 'PR'].set_index('Plan Name') + +# Melt the original and PR dataframes +original_melted = original_df.reset_index().melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='Original') +pr_melted = pr_df.reset_index().melt(id_vars=['Plan Name'], value_vars=['3%3mm', '2%3mm', '1%3mm'], var_name='Metric', value_name='PR') + +metrics = ['3%3mm', '2%3mm', '1%3mm'] + +for i, metric in enumerate(metrics): + # Create the plot for each metric + fig, ax = plt.subplots(figsize=(10, 8)) + sns.set(style="whitegrid") + + # Plot the original values + subset_df = original_melted[original_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['Original'], y=subset_df['Plan Name'], color='darkgrey', legend=False, s=100, marker='o', alpha=0.5, ax=ax) + + # Plot the PR values + subset_df = pr_melted[pr_melted['Metric'] == metric] + sns.scatterplot(x=subset_df['PR'], y=subset_df['Plan Name'], color='darkblue', legend=False, s=100, marker='o', alpha=0.5, ax=ax) + + # Set the limits for the x-axis and set label color to dark grey + ax.set_xlim([50, 110]) + ax.set_xlabel('Pass Rates', color='darkgrey') + + # Add labels and title + ax.set_ylabel('Plan Name') + ax.set_title(f'PR {metric}', loc='left') # Set title to just the metric and align left + + plt.tight_layout() + plt.savefig(f'PR_{metric}_pass_rate_plot.png', dpi=300) + plt.show() \ No newline at end of file diff --git a/src/plots_for_presentation/legend.py b/src/plots_for_presentation/legend.py new file mode 100644 index 0000000..509ab7f --- /dev/null +++ b/src/plots_for_presentation/legend.py @@ -0,0 +1,19 @@ +import matplotlib.pyplot as plt +from matplotlib.lines import Line2D + +# Create legend handles manually +legend_elements = [Line2D([0], [0], marker='o', color='w', markerfacecolor='darkblue', markersize=10, label='Pass Rate Corrected'), + Line2D([0], [0], marker='o', color='w', markerfacecolor='darkgrey', markersize=10, label='Pass Rate Uncorrected'), + Line2D([0], [0], marker='^', color='w', markerfacecolor='green', markersize=10, label='Pass Rate Difference Better'), + Line2D([0], [0], marker='^', color='w', markerfacecolor='red', markersize=10, label='Pass Rate Difference Worse')] + +# Create the figure and the axes +fig, ax = plt.subplots() + +# Add the legend to the plot +ax.legend(handles=legend_elements, loc='center') + +# Hide axes +ax.axis('off') +plt.savefig(f'legend.png', dpi=300) +plt.show() \ No newline at end of file From c793880acf4977ce86574a15a091c703e82e06be Mon Sep 17 00:00:00 2001 From: kkaan Date: Mon, 17 Jun 2024 13:32:36 +1000 Subject: [PATCH 7/8] added a legend --- src/plots_for_presentation/difference_plot_dpp.py | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/src/plots_for_presentation/difference_plot_dpp.py b/src/plots_for_presentation/difference_plot_dpp.py index e2bb8b2..68121ac 100644 --- a/src/plots_for_presentation/difference_plot_dpp.py +++ b/src/plots_for_presentation/difference_plot_dpp.py @@ -41,6 +41,13 @@ ax.set_ylabel('Plan Name') ax.set_title(f'DPP {metric}', loc='left') # Set title to just the metric and align left + # Add legend + legend_labels = ['Original', 'DPP'] + legend_colors = ['darkgrey', 'darkblue'] + patches = [plt.plot([],[], marker="o", ms=10, ls="", mec=None, color=legend_colors[i], + label="{:s}".format(legend_labels[i]) )[0] for i in range(len(legend_labels))] + ax.legend(handles=patches, bbox_to_anchor=(0, 0.5), loc='center left', ncol=1) + plt.tight_layout() plt.savefig(f'DPP_{metric}_pass_rate_plot.png', dpi=300) plt.show() \ No newline at end of file From e5dd1b67a120e2d12cf805645d2d9267bdf42b52 Mon Sep 17 00:00:00 2001 From: kkaan Date: Thu, 20 Jun 2024 14:44:16 +1000 Subject: [PATCH 8/8] added a legend --- src/plots_for_presentation/difference_plot_dpp.py | 5 +++-- src/plots_for_presentation/difference_plot_pr.py | 10 +++++++++- 2 files changed, 12 insertions(+), 3 deletions(-) diff --git a/src/plots_for_presentation/difference_plot_dpp.py b/src/plots_for_presentation/difference_plot_dpp.py index 68121ac..6b9e8b7 100644 --- a/src/plots_for_presentation/difference_plot_dpp.py +++ b/src/plots_for_presentation/difference_plot_dpp.py @@ -24,7 +24,8 @@ # Create the plot for each metric fig, ax = plt.subplots(figsize=(10, 8)) sns.set(style="whitegrid") - + # Ensure grid is turned on + ax.grid(True) # Plot the original values subset_df = original_melted[original_melted['Metric'] == metric] sns.scatterplot(x=subset_df['Original'], y=subset_df['Plan Name'], color='darkgrey', legend=False, s=100, marker='o', alpha=0.5, ax=ax) @@ -42,7 +43,7 @@ ax.set_title(f'DPP {metric}', loc='left') # Set title to just the metric and align left # Add legend - legend_labels = ['Original', 'DPP'] + legend_labels = ['Original', 'DPP Corrected'] legend_colors = ['darkgrey', 'darkblue'] patches = [plt.plot([],[], marker="o", ms=10, ls="", mec=None, color=legend_colors[i], label="{:s}".format(legend_labels[i]) )[0] for i in range(len(legend_labels))] diff --git a/src/plots_for_presentation/difference_plot_pr.py b/src/plots_for_presentation/difference_plot_pr.py index 0ab42e7..fe5168a 100644 --- a/src/plots_for_presentation/difference_plot_pr.py +++ b/src/plots_for_presentation/difference_plot_pr.py @@ -31,7 +31,7 @@ # Plot the PR values subset_df = pr_melted[pr_melted['Metric'] == metric] - sns.scatterplot(x=subset_df['PR'], y=subset_df['Plan Name'], color='darkblue', legend=False, s=100, marker='o', alpha=0.5, ax=ax) + sns.scatterplot(x=subset_df['PR'], y=subset_df['Plan Name'], color='green', legend=False, s=100, marker='o', alpha=0.5, ax=ax) # Set the limits for the x-axis and set label color to dark grey ax.set_xlim([50, 110]) @@ -41,6 +41,14 @@ ax.set_ylabel('Plan Name') ax.set_title(f'PR {metric}', loc='left') # Set title to just the metric and align left + # Add legend + legend_labels = ['Original', 'PR Corrected'] + legend_colors = ['darkgrey', 'green'] + patches = [plt.plot([],[], marker="o", ms=10, ls="", mec=None, color=legend_colors[i], + label="{:s}".format(legend_labels[i]) )[0] for i in range(len(legend_labels))] + ax.legend(handles=patches, bbox_to_anchor=(0, 0.5), loc='center left', ncol=1) + + plt.tight_layout() plt.savefig(f'PR_{metric}_pass_rate_plot.png', dpi=300) plt.show() \ No newline at end of file