Splitting reading and interpolating into 2 functions

Change-Id: I02ada28a6a04cc6dd91dae1e560eb72ceb3f8987
This commit is contained in:
t2 2017-10-10 12:20:15 +02:00 committed by Tomasz Dyczek
parent d0f71e336c
commit aa25d38639
7 changed files with 101 additions and 64 deletions

View File

@ -2,12 +2,17 @@ import matplotlib.pyplot as plt
import csv import csv
import numpy as np import numpy as np
from numpy import linalg from numpy import linalg
from scipy.interpolate import interp1d
class NoInterpolatedDataException(Exception): class NoInterpolatedDataException(Exception):
pass pass
class WrongWavelengthsInConf(Exception):
pass
class OpusData(object): class OpusData(object):
def __init__(self, raw_data, interp=None, meta=None): def __init__(self, raw_data, interp=None, meta=None):
self.raw_data = raw_data self.raw_data = raw_data
@ -46,3 +51,39 @@ class OpusData(object):
self.interpolated_data[1], self.interpolated_data[1],
'ro') 'ro')
plt.show() 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

View File

@ -1,8 +1,7 @@
import numpy as np import numpy as np
from struct import unpack_from from struct import unpack_from
from scipy.interpolate import interp1d from .utils import find_all
from utils import find_all from .opus_data import OpusData
from opus_data import OpusData
from itertools import product from itertools import product
@ -10,11 +9,7 @@ class NoAbsorbanceSpectra(Exception):
pass pass
class WrongWavelengthsInConf(Exception): def opus_reader(filepath):
pass
def opus_reader(filepath, conf, interpolate=True):
with open(filepath, 'rb') as f: with open(filepath, 'rb') as f:
buff = f.read() buff = f.read()
@ -23,19 +18,11 @@ def opus_reader(filepath, conf, interpolate=True):
# Choose best ab spectra # Choose best ab spectra
ab_spectra, ab_wavenumbers = choose_ab(fxv_spc, spc, wavenumbers) 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) meta = get_meta_data(buff)
if interpolate: return OpusData(zip(*wave_num_abs_pair), meta=meta)
return OpusData(zip(*pairs), interp=(iwavenumber, yi), meta=meta)
else:
return OpusData(zip(*pairs), meta=meta)
def choose_ab(fxv_spc, spc, wavenumbers): def choose_ab(fxv_spc, spc, wavenumbers):
@ -123,44 +110,6 @@ def generate_wavelengths(lxv_spc, fxv_spc, npt_spc):
return wavenumbers 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): def get_meta_data(buff):
# Getting source of instruments # Getting source of instruments
all_ins = tuple(find_all('INS', buff)) all_ins = tuple(find_all('INS', buff))

26
example.py Normal file
View File

@ -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])

View File

@ -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=""))
}

3
requirements.in Normal file
View File

@ -0,0 +1,3 @@
numpy
scipy
matplotlib

View File

@ -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

10
setup.py Normal file
View File

@ -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'],
)