diff --git a/diffpy/labpdfproc/functions.py b/diffpy/labpdfproc/functions.py index 2db6da3..b111028 100644 --- a/diffpy/labpdfproc/functions.py +++ b/diffpy/labpdfproc/functions.py @@ -4,13 +4,32 @@ def compute_cve(diffraction_data, mud, wavelength): # for a given mu and d and lambda, we will compute cve on a tth grid # something arbitrary for the moment - cve = Diffraction_object() + + # calculate corresponding cve for the given mud + radius_mm = 1 + n_points_on_diameter = 249 + mu = mud/2 + tth_grid = np.arange(1, 141, 1) + abs_correction = Gridded_circle(radius_mm, n_points_on_diameter, mu=mu) + distances, muls = [], [] + + for angle in tth_grid: + abs_correction.set_distances_at_angle(angle) + abs_correction.set_muls_at_angle(angle) + distances.append(sum(abs_correction.distances)) + muls.append(sum(abs_correction.muls)) + distances = np.array(distances) / abs_correction.total_points_in_grid + muls = np.array(muls) / abs_correction.total_points_in_grid + + # interpolate + cve = Diffraction_object(wavelength=wavelength) cve_x = diffraction_data.on_tth[0] - cve_y = cve_x * mud * wavelength - cve.insert_scattering_quantity(cve_x, cve_y, "tth", metadata={ }) + cve_y = np.interp(cve_x, tth_grid, muls) + cve.insert_scattering_quantity(cve_x, cve_y, "tth", metadata={}) return cve + def apply_corr(i_m, cve): # we apply the absorption correction by doing: I(tth) * c_ve - i_c = i_m * cve + i_c = i_m / cve.on_tth[1] return i_c diff --git a/diffpy/labpdfproc/labpdfprocapp.py b/diffpy/labpdfproc/labpdfprocapp.py index 7977b58..78f45f0 100644 --- a/diffpy/labpdfproc/labpdfprocapp.py +++ b/diffpy/labpdfproc/labpdfprocapp.py @@ -1,38 +1,16 @@ import sys from argparse import ArgumentParser +import numpy as np from diffpy.labpdfproc.functions import compute_cve, apply_corr, wavelengths from diffpy.utils.parsers.loaddata import loadData from diffpy.utils.scattering_objects.diffraction_objects import Diffraction_object known_sources = ["Ag", "Mo"] -# def load_data(input_file): -# # we want to load .xy, xye, file types. These are the most common. For anyting else (.snc, .txt, .csv, .dat) we return an error message. Why: some of them have different delimineters. -# # we want to load the tth column in a numpy array -# # we want to load the intensitie Im into a numpy array -# # the input files should not contain any header. Typically, .xy or .xye files don't contain headers. -# tth = np.loadtxt(input_file, usecols=0) -# i_m = np.loadtxt(input_file, usecols=1) -# # this should return an error if the first row contains anything except a float, and if the columns are not separated by a space. -# # I think the latter is also dealt with if we check if the first elemnt in one tth is a flaot. -# if np.issubdtype(tth[0], np.floating) and np.issubdtype(i_m[0], np.floating): -# return tth, i_m -# else: -# raise ValueError('Error: your .xy contains headers. Delete the header rows in your .xy or .xye file') -# def tth_to_q(tth, wl): -# tth_rad = np.deg2rad(tth) -# q = (4 * np.pi / wl) * np.sin(tth_rad / 2) -# return q - -# def write_files(base_name): - # we save the corrected intensities in a two-column file on a q-grid and a tth-grid. - # we need to know the x-ray wavelenth so that we can convert tth to q - # we make a new two-column file .chi where column 1 contains the q grid and columnt 2 contains the corrected intensities. - def get_args(): p = ArgumentParser() p.add_argument("data_file", help="the filename of the datafile to load") @@ -51,20 +29,15 @@ def main(): args = get_args() wavelength = wavelengths[args.anode_type] - input_pattern = Diffraction_object() + input_pattern = Diffraction_object(wavelength=wavelength) xarray, yarray = loadData(args.data_file, unpack=True) input_pattern.insert_scattering_quantity(xarray, yarray, "tth", metadata={ }) cve = compute_cve(input_pattern, args.mud, wavelength) - i_c = input_pattern * cve + i_c = apply_corr(yarray, cve) - # get the basename from the input_file and save the corrected patter as a .tth and a .chi file. - # base_name = input_file.split('.')[0] - # output_chi = f"{base_name}.chi" - # output_tth = f"{base_name}.tth" - # np.savetxt(output_tth, np.column_stack((tth, i_c)), header='tth I(tth)') - # np.savetxt(output_chi, np.column_stack((q, i_c)), header='tth I(tth)') - input_pattern.dump("filename", type="chi") + data_to_save = np.column_stack((xarray, i_c)) + np.savetxt(f'{base_name}_proc.chi', data_to_save) if __name__ == '__main__': main()