Commit 31968a5a authored by Zhang, Yuanpeng's avatar Zhang, Yuanpeng
Browse files

change peak shape to lorentzian

parent 63784a82
Loading
Loading
Loading
Loading
+20697 −20697

File changed.

Preview size limit exceeded, changes collapsed.

+4 −4
Original line number Diff line number Diff line
{
    "Diamond": "/SNS/NOM/IPTS-28922/nexus/NOM_172394.nxs.h5",
    "Diamond": "/SNS/NOM/IPTS-30147/nexus/NOM_189201.nxs.h5",
    "Instrument": "NOM",
    "Date": "2022-05-25",
    "SampleEnv": "shifter",
    "Date": "2023-01-10",
    "SampleEnv": "furnace",
    "OutputDir": "/SNS/NOM/shared/autoreduce/calibration",
    "GenShadowMask": "shadow_mask_shifter_172394.in",
    "GenShadowMask": "shadow_mask_shifter_189201.in",
    "DiaLattParam": 3.5671299351,
    "GroupMethod": "KMEANS_ED",
    "SaveInitCalTable": true,
+109 −13
Original line number Diff line number Diff line
@@ -570,7 +570,7 @@ for bank in range(mtd["dia_calibrant_d_cc_f_tof"].getNumberHistograms()):
        if bank == 0:
            chi_s_min = 5.E-3
        else:
            chi_s_min = 5.E-4
            chi_s_min = 5.E-3

        d_norm = list()
        d_fitted = list()
@@ -590,9 +590,9 @@ for bank in range(mtd["dia_calibrant_d_cc_f_tof"].getNumberHistograms()):
                    OutputWorkspace="dia_calibrant_d_cc_f_fine_fit" + str(bank + 1) + str(key),
                    ParameterTableWorkspace="dia_calibrant_d_cc_f_fine_table" + str(bank + 1) + str(key),
                    WorkspaceIndex=bank,
                    PeakFunctionType="BackToBackExponential",
                    PeakParameterNames="I, A, B, X0, S",
                    PeakParameterValues=f"{height_est}, 100, 100, {c_est}, 0.0005",
                    PeakFunctionType="Lorentzian",
                    PeakParameterNames="Amplitude, PeakCentre, FWHM",
                    PeakParameterValues=f"{height_est}, {c_est}, 0.05",
                    BackgroundType="Linear",
                    BackgroundParameterNames="A0, A1",
                    BackgroundParameterValues="10, 0",
@@ -653,7 +653,7 @@ for bank in range(mtd["dia_calibrant_d_cc_f_tof"].getNumberHistograms()):

        PDCalibration(InputWorkspace=f"dia_calibrant_d_cc_f_tof_{bank + 1}",
                      TofBinning="300, -0.0001, 16667",
                      PeakFunction="BackToBackExponential",
                      PeakFunction="Lorentzian",
                      BackgroundType="Linear",
                      PeakPositions=center_str,
                      PeakWindow=range_str,
@@ -702,9 +702,9 @@ for bank in range(mtd["dia_calibrant_d_cc_f_tof"].getNumberHistograms()):
                                OutputWorkspace=f"dia_calibrant_d_cc_f_tof_{bank + 1}_{row_i}_fit",
                                ParameterTableWorkspace=f"dia_calibrant_d_cc_f_tof_{bank + 1}_{row_i}_table",
                                WorkspaceIndex=0,
                                PeakFunctionType="BackToBackExponential",
                                PeakParameterNames="I, A, B, X0, S",
                                PeakParameterValues=f"{height_est}, 100., 100., {tof_cent_trial}, 0.0005",
                                PeakFunctionType="Lorentzian",
                                PeakParameterNames="Amplitude, PeakCentre, FWHM",
                                PeakParameterValues=f"{height_est}, {tof_cent_trial}, 100",
                                BackgroundType="Linear",
                                BackgroundParameterNames="A0, A1",
                                BackgroundParameterValues="10, 0",
@@ -730,9 +730,9 @@ for bank in range(mtd["dia_calibrant_d_cc_f_tof"].getNumberHistograms()):
                            OutputWorkspace=f"dia_calibrant_d_cc_f_tof_{bank + 1}_{row_i}_fit",
                            ParameterTableWorkspace=f"dia_calibrant_d_cc_f_tof_{bank + 1}_{row_i}_table",
                            WorkspaceIndex=0,
                            PeakFunctionType="BackToBackExponential",
                            PeakParameterNames="I, A, B, X0, S",
                            PeakParameterValues=f"{height_est}, 100., 100., {tof_cent_opt}, 0.0005",
                            PeakFunctionType="Lorentzian",
                            PeakParameterNames="Amplitude, PeakCentre, FWHM",
                            PeakParameterValues=f"{height_est}, {tof_cent_opt}, 100",
                            BackgroundType="Linear",
                            BackgroundParameterNames="A0, A1",
                            BackgroundParameterValues="10, 0",
@@ -772,6 +772,14 @@ for bank in range(mtd["dia_calibrant_d_cc_f_tof"].getNumberHistograms()):
            k = float(key.split("-")[1])
            l = float(key.split("-")[2])
            ds_tmp = dia_a / np.sqrt(h**2. + k**2. + l**2.)
            if bank == 3:
                if ds_tmp > 1.1:
                    break
            elif bank == 4:
                if ds_tmp > 0.7:
                    break
            else:
                pass
            center_str.append("{0:8.6F}".format(ds_tmp))
            range_str.append("{0:8.6F}".format(item["dmin"]))
            range_str.append("{0:8.6F}".format(item["dmax"]))
@@ -789,9 +797,9 @@ for bank in range(mtd["dia_calibrant_d_cc_f_tof"].getNumberHistograms()):

        if bank == 5:
            range_str = "0.3"
            pf_use = "Gaussian"
            pf_use = "Lorentzian"
        else:
            pf_use = "BackToBackExponential"
            pf_use = "Lorentzian"

        PDCalibration(InputWorkspace=f"dia_calibrant_d_cc_f_tof_{bank + 1}",
                      TofBinning="300, -0.0001, 16667",
@@ -845,6 +853,94 @@ for bank in range(mtd["dia_calibrant_d_cc_f_tof"].getNumberHistograms()):
                d_tmp /= difc_arb
                d_fitted.append(d_tmp)

        if bank in (3, 4):
            center_str = list()
            range_str = list()
            row_count = 0
            peak_list = dict()
            if bank == 4:
                chi_s_min = 5.E-4
                chi_s_min_pd = 5.E-5
            else:
                chi_s_min = 5.E-3
                chi_s_min_pd = 5.E-3

            for key, item in dia_params[f"{bank + 1}"].items():
                h = float(key.split("-")[0])
                k = float(key.split("-")[1])
                l = float(key.split("-")[2])
                ds_tmp = dia_a / np.sqrt(h**2. + k**2. + l**2.)
                if bank == 3:
                    if ds_tmp < 1.1:
                        continue
                elif bank == 4:
                    if ds_tmp < 0.7:
                        continue
                else:
                    pass
                center_str.append("{0:8.6F}".format(ds_tmp))
                range_str.append("{0:8.6F}".format(item["dmin"]))
                range_str.append("{0:8.6F}".format(item["dmax"]))
                peak_list[row_count] = ["{0:8.6F}".format(ds_tmp),
                                        "{0:8.6F}".format(item["dmin"]),
                                        "{0:8.6F}".format(item["dmax"])]
                row_count += 1

            center_str = ",".join(center_str)
            range_str = ",".join(range_str)

            pf_use = "Lorentzian"

            PDCalibration(InputWorkspace=f"dia_calibrant_d_cc_f_tof_{bank + 1}",
                          TofBinning="300, -0.0001, 16667",
                          PeakFunction=pf_use,
                          BackgroundType="Linear",
                          PeakPositions=center_str,
                          PeakWindow=range_str,
                          CalibrationParameters="DIFC+TZERO+DIFA",
                          OutputCalibrationTable=f"pdcal_table_{bank + 1}_high_d",
                          DiagnosticWorkspaces=f"pdcal_diag_{bank + 1}_high_d")

            difc_arb = float(mtd["_det_f_table"].row(bank)["DIFC"])
            for row_i in range(mtd[f"pdcal_diag_{bank + 1}_high_d_fitparam"].rowCount()):
                chi2_tmp = mtd[f"pdcal_diag_{bank + 1}_high_d_fitparam"].row(row_i)["chi2"]
                int_tmp = mtd[f"pdcal_diag_{bank + 1}_high_d_fitparam"].row(row_i)["intensity"]
                if row_i > 0:
                    w_current = mtd[f"pdcal_diag_{bank + 1}_high_d_fitparam"].row(row_i)["width"]
                    w_last = mtd[f"pdcal_diag_{bank + 1}_high_d_fitparam"].row(row_i - 1)["width"]
                    h_current = mtd[f"pdcal_diag_{bank + 1}_high_d_fitparam"].row(row_i)["height"]
                    h_last = mtd[f"pdcal_diag_{bank + 1}_high_d_fitparam"].row(row_i - 1)["height"]
                    extra_cond = abs(w_last - w_current) < 1E-10 and abs(h_current - h_last) < 1E-10
                else:
                    extra_cond = int_tmp < 0 or chi2_tmp / int_tmp > chi_s_min_pd

                if int_tmp < 0 or chi2_tmp / int_tmp > chi_s_min_pd or extra_cond:
                    second_peak_right = "{0:8.6F}".format(
                        2. * float(peak_list[row_i][2]) - float(peak_list[row_i][0])
                    )
                    PDCalibration(InputWorkspace=f"dia_calibrant_d_cc_f_tof_{bank + 1}",
                          TofBinning="300, -0.0001, 16667",
                          PeakFunction=pf_use,
                          BackgroundType="Linear",
                          PeakPositions=(float(peak_list[row_i][0]), float(peak_list[row_i][2])),
                          PeakWindow=f"{peak_list[row_i][1]}, {peak_list[row_i][2]}, {peak_list[row_i][0]}, {second_peak_right}",
                          CalibrationParameters="DIFC+TZERO+DIFA",
                          OutputCalibrationTable=f"pdcal_table_{bank + 1}_{row_i}_high_d",
                          DiagnosticWorkspaces=f"pdcal_diag_{bank + 1}_{row_i}_high_d")

                    chi2_tmp = mtd[f"pdcal_diag_{bank + 1}_{row_i}_high_d_fitparam"].row(0)["chi2"]
                    int_tmp = mtd[f"pdcal_diag_{bank + 1}_{row_i}_high_d_fitparam"].row(0)["intensity"]
                    if int_tmp > 0 and chi2_tmp / int_tmp < chi_s_min:
                        d_norm.append(float(peak_list[row_i][0]))
                        d_tmp = mtd[f"pdcal_diag_{bank + 1}_{row_i}_high_d_fitparam"].row(0)["centre"]
                        d_tmp /= difc_arb
                        d_fitted.append(d_tmp)
                else:
                    d_norm.append(float(peak_list[row_i][0]))
                    d_tmp = mtd[f"pdcal_diag_{bank + 1}_high_d_fitparam"].row(row_i)["centre"]
                    d_tmp /= difc_arb
                    d_fitted.append(d_tmp)

        d_norm.sort()
        d_fitted.sort()

utils/nom_cal_btb.py

0 → 100644
+1179 −0

File added.

Preview size limit exceeded, changes collapsed.

+1275 −0

File added.

Preview size limit exceeded, changes collapsed.

Loading