Commit fb53f848 authored by Zhang, Yuanpeng's avatar Zhang, Yuanpeng ⭐
Browse files

add in alignment diagnostics

parent 032626b7

utils/pg3_calib_doctor.py

+54 −0
Original line number Diff line number Diff line
@@ -208,6 +208,18 @@ def pg3_calib_doctor(run_num, cal_file, manual_mask, out_dir, diag_dict): 
        OutputWorkspace="manual_mask"

    )



    pg3_chunks = [

        [[0, 15000], [0.47, 0.49], [0.4768, 0]],

        [[15000, 19580], [0.62, 0.65], [0.6307, 0]],

        [[19580, 25000], [0.66, 0.71], [0.6867, 0]],

        [[25000, 29057], [1.00, 1.15], [1.0758, 0]],

        [[29057, 30184], [1.80, 2.30], [2.0600, 0]],

        [[30184, 33431], [0.47, 0.49], [0.4768, 0]],

        [[33431, 35567], [0.55, 0.58], [0.5642, 0]],

        [[35567, 39000], [0.66, 0.71], [0.6867, 0]],

        [[39000, 43120], [1.25, 1.275], [1.2615, 0]]

    ]



    pg3_chunks_1 = [

        [[0, 9750], [0.4, 0.55], 9],

        [[9750, 25000], [0.5, 1.0], 9],
@@ -244,6 +256,7 @@ def pg3_calib_doctor(run_num, cal_file, manual_mask, out_dir, diag_dict): 
    bad_pixels_unmasked = []

    num_hist = mtd["pg3_input"].getNumberHistograms()

    num_masks = 0

    rel_offset = []

    for i in range(num_hist):

        idx = get_chunk_index(i)

        d_range_1 = pg3_chunks_1[idx][1]
@@ -262,6 +275,26 @@ def pg3_calib_doctor(run_num, cal_file, manual_mask, out_dir, diag_dict): 
        if mtd["_mask"].readY(i)[0] == 1.0:

            num_masks += 1



        # Alignment diagnostics

        if mtd["_mask"].readY(i)[0] == 1.0:

            rel_offset.append(0.0)

            continue



        x = mtd["pg3_wksp_d"].readX(i)

        y = mtd["pg3_wksp_d"].readY(i)

        y = np.insert(y, 0, y[0])



        x_min, x_max = pg3_chunks[idx][1][0], pg3_chunks[idx][1][1]



        mask = (x >= x_min) & (x <= x_max)

        x_peak = x[mask]

        y_peak = y[mask]



        y_peak = np.where(y_peak < 0.03 * y_peak.max(), 0, y_peak)

        center_of_mass = np.sum(x_peak * y_peak) / np.sum(y_peak)



        rel_offset.append((center_of_mass - pg3_chunks[idx][2][0]) / center_of_mass)



    print("[Info] # of masked good pixels: ", len(good_pixels_masked))

    print("[Info] # of unmasked bad pxiels: ", len(bad_pixels_unmasked))
@@ -298,6 +331,27 @@ def pg3_calib_doctor(run_num, cal_file, manual_mask, out_dir, diag_dict): 
        diag_input = [num_masks, len(good_pixels_masked)]

        plot_masked_spectra(run_num, run_out_dir, False, diag_dict, diag_input)



    data = np.array(rel_offset)

    data = np.abs(data)

    original_indices = np.where(data != 0)[0]

    data = data[original_indices]



    # --- Outlier detection using IQR method ---

    Q1 = np.percentile(data, 5)

    Q3 = np.percentile(data, 95)

    IQR = Q3 - Q1

    lower = Q1 - 1.5 * IQR

    upper = Q3 + 1.5 * IQR



    mean = np.mean(data)

    std = np.std(data)



    outlier_mask = (data < lower) | (data > upper)

    outliers = data[outlier_mask]

    outlier_indices_original = original_indices[outlier_mask]



    return





_DIAGNOSTICS_KEYS = ["JANISGAS", "LTJANIS", "MICAS", "PAC", "MAG", "OC"]