From aa25d3863988ac3b195d7437f5d408e2087b1058 Mon Sep 17 00:00:00 2001 From: t2 Date: Tue, 10 Oct 2017 12:20:15 +0200 Subject: [PATCH] Splitting reading and interpolating into 2 functions Change-Id: I02ada28a6a04cc6dd91dae1e560eb72ceb3f8987 --- brukeropusreader/opus_data.py | 41 ++++++++++++++++++++++ brukeropusreader/opus_reader.py | 61 +++------------------------------ example.py | 26 ++++++++++++++ get_bart_result.R | 6 ---- requirements.in | 3 ++ requirements.txt | 18 ++++++++-- setup.py | 10 ++++++ 7 files changed, 101 insertions(+), 64 deletions(-) create mode 100644 example.py delete mode 100644 get_bart_result.R create mode 100644 requirements.in create mode 100644 setup.py diff --git a/brukeropusreader/opus_data.py b/brukeropusreader/opus_data.py index 3e8decc..7c4ab27 100644 --- a/brukeropusreader/opus_data.py +++ b/brukeropusreader/opus_data.py @@ -2,12 +2,17 @@ import matplotlib.pyplot as plt import csv import numpy as np from numpy import linalg +from scipy.interpolate import interp1d class NoInterpolatedDataException(Exception): pass +class WrongWavelengthsInConf(Exception): + pass + + class OpusData(object): def __init__(self, raw_data, interp=None, meta=None): self.raw_data = raw_data @@ -46,3 +51,39 @@ class OpusData(object): self.interpolated_data[1], 'ro') plt.show() + + def gen_interp_waves(self, conf): + wunit = (conf.wave_start - conf.wave_end) / (float(conf.iwsize - 1)) + a = conf.wave_start + wunit + iwavenumber = [0] * conf.iwsize + for i in xrange(len(iwavenumber)): + a -= wunit + iwavenumber[i] = a + return iwavenumber + + def interpolate_spectra(self, conf): + xav, yav = self.raw_data[0], self.raw_data[1] + ab_wavenumbers = self.raw_data[0] + iwavenumber = self.gen_interp_waves(conf) + + xa_min = min(ab_wavenumbers) + xa_max = max(ab_wavenumbers) + n_interp = conf.iwsize + yi = [] + + f2 = interp1d(xav, yav, kind='cubic', assume_sorted=True) + + for k in range(n_interp): + if iwavenumber[k] < xa_min: + if k == n_interp - 1: + yi.append(yav[0]) + else: + print("Wrong wavelenghts for interpolating data") + elif iwavenumber[k] > xa_max: + if k == 0: + yi.append(yav[-1]) + else: + raise WrongWavelengthsInConf() + else: + yi.append(f2(iwavenumber[k])) + return iwavenumber, yi diff --git a/brukeropusreader/opus_reader.py b/brukeropusreader/opus_reader.py index 096d754..bf09b10 100644 --- a/brukeropusreader/opus_reader.py +++ b/brukeropusreader/opus_reader.py @@ -1,8 +1,7 @@ import numpy as np from struct import unpack_from -from scipy.interpolate import interp1d -from utils import find_all -from opus_data import OpusData +from .utils import find_all +from .opus_data import OpusData from itertools import product @@ -10,11 +9,7 @@ class NoAbsorbanceSpectra(Exception): pass -class WrongWavelengthsInConf(Exception): - pass - - -def opus_reader(filepath, conf, interpolate=True): +def opus_reader(filepath): with open(filepath, 'rb') as f: buff = f.read() @@ -23,19 +18,11 @@ def opus_reader(filepath, conf, interpolate=True): # Choose best ab spectra ab_spectra, ab_wavenumbers = choose_ab(fxv_spc, spc, wavenumbers) - pairs = reversed(zip(ab_wavenumbers, ab_spectra)) + wave_num_abs_pair = reversed(zip(ab_wavenumbers, ab_spectra)) - # AB spectra and wave extracted, now interpolating - if interpolate: - iwavenumber, yi = interpolate_spectra(pairs, ab_wavenumbers, conf) - - # Spectra succesfully interpolated, now extracting metadata meta = get_meta_data(buff) - if interpolate: - return OpusData(zip(*pairs), interp=(iwavenumber, yi), meta=meta) - else: - return OpusData(zip(*pairs), meta=meta) + return OpusData(zip(*wave_num_abs_pair), meta=meta) def choose_ab(fxv_spc, spc, wavenumbers): @@ -123,44 +110,6 @@ def generate_wavelengths(lxv_spc, fxv_spc, npt_spc): return wavenumbers -def generate_interpolated_wavenumbers(conf): - wunit = (conf.wave_start - conf.wave_end) / (float(conf.iwsize - 1)) - a = conf.wave_start + wunit - iwavenumber = [0] * conf.iwsize - for i in xrange(len(iwavenumber)): - a -= wunit - iwavenumber[i] = a - return iwavenumber - - -def interpolate_spectra(pairs, ab_wavenumbers, conf): - iwavenumber = generate_interpolated_wavenumbers(conf) - - xav, yav = zip(*pairs)[0], zip(*pairs)[1] - - xa_min = min(ab_wavenumbers) - xa_max = max(ab_wavenumbers) - n_interp = conf.iwsize - yi = [] - - f2 = interp1d(xav, yav, kind='cubic', assume_sorted=True) - - for k in range(n_interp): - if iwavenumber[k] < xa_min: - if k == n_interp - 1: - yi.append(yav[0]) - else: - print("Wrong wavelenghts for interpolating data") - elif iwavenumber[k] > xa_max: - if k == 0: - yi.append(yav[-1]) - else: - raise WrongWavelengthsInConf() - else: - yi.append(f2(iwavenumber[k])) - return iwavenumber, yi - - def get_meta_data(buff): # Getting source of instruments all_ins = tuple(find_all('INS', buff)) diff --git a/example.py b/example.py new file mode 100644 index 0000000..014f3ee --- /dev/null +++ b/example.py @@ -0,0 +1,26 @@ +from brukeropusreader.opus_reader import opus_reader +from os import walk +import os +import sys + + +def filter_fun(nm): + return (nm[0] != '.') and (nm[-2:] in ('.0', '.1', '.2', '.3')) + + +def main(path_to_opus_files): + i = 0 + all_files = 0 + for path, subdirs, files in walk(path_to_opus_files): + for name in filter(filter_fun, files): + all_files += 1 + try: + od = opus_reader(os.path.join(path, name)) + od.plot_raw() + except Exception as e: + i += 1 + print i, ":", name, str(e) + + +if __name__ == '__main__': + main(sys.argv[1]) diff --git a/get_bart_result.R b/get_bart_result.R deleted file mode 100644 index cb435d5..0000000 --- a/get_bart_result.R +++ /dev/null @@ -1,6 +0,0 @@ -library(rbart) - -save.soil <- function(opus.file.path, conf.path, output.file='result.csv', output.path='./'){ - df = read.opus(opus.file.path, conf.path) - write.csv(df, paste(output.path, output.file, sep="")) -} \ No newline at end of file diff --git a/requirements.in b/requirements.in new file mode 100644 index 0000000..f23634b --- /dev/null +++ b/requirements.in @@ -0,0 +1,3 @@ +numpy +scipy +matplotlib diff --git a/requirements.txt b/requirements.txt index 6bad103..f69a35d 100644 --- a/requirements.txt +++ b/requirements.txt @@ -1,2 +1,16 @@ -numpy -scipy +# +# This file is autogenerated by pip-compile +# To update, run: +# +# pip-compile --output-file requirements.txt requirements.in +# +backports.functools-lru-cache==1.4 # via matplotlib +cycler==0.10.0 # via matplotlib +matplotlib==2.1.0 +numpy==1.13.3 +pyparsing==2.2.0 # via matplotlib +python-dateutil==2.6.1 # via matplotlib +pytz==2017.2 # via matplotlib +scipy==0.19.1 +six==1.11.0 # via cycler, matplotlib, python-dateutil +subprocess32==3.2.7 # via matplotlib diff --git a/setup.py b/setup.py new file mode 100644 index 0000000..896d6bd --- /dev/null +++ b/setup.py @@ -0,0 +1,10 @@ +from distutils.core import setup + + +setup(name='bruker-opus-reader', + version='1.0', + description='Bruker OPUS files reader', + author='t2', + author_email='t2@qed.ai', + packages=['brukeropusreader'], + )