What a nice try! Let's fight the pytest!

TESTTTTTT/new.py

@@ -0,0 +1,3 @@
+from src import *
+print(mod1.f_mod1())
+print(mod2.f_mod2())

TESTTTTTT/src/mod1.py

@@ -0,0 +1,2 @@
+def f_mod1():
+    return 1

TESTTTTTT/src/mod2.py

@@ -0,0 +1,3 @@
+from mod1 import f_mod1
+def f_mod2():
+    return f_mod1()

TESTTTTTT/tests/test.py

@@ -0,0 +1,8 @@
+from new.mod1 import f_mod1
+from new.mod2 import f_mod2
+
+def test_f_mod1():
+    assert f_mod1() == 1
+
+def test_f_mod2():
+    assert f_mod2() == 1

requirements.txt

@@ -0,0 +1,86 @@
+apturl==0.5.2
+attrs==22.2.0
+BaselineRemoval==0.1.1
+blinker==1.4
+brukeropusreader==1.3.4
+catfish==4.16.3
+cffi==1.15.0
+chardet==4.0.0
+click==8.0.3
+colorama==0.4.4
+command-not-found==0.3
+contourpy==1.0.5
+coverage==7.1.0
+cryptography==3.4.8
+cupshelpers==1.0
+cycler==0.11.0
+dbus-python==1.2.18
+defer==1.0.6
+distro==1.7.0
+distro-info===1.1build1
+exceptiongroup==1.1.0
+fonttools==4.37.4
+httplib2==0.20.2
+importlib-metadata==4.6.4
+iniconfig==2.0.0
+jeepney==0.7.1
+joblib==1.2.0
+keyring==23.5.0
+kiwisolver==1.4.4
+language-selector==0.1
+launchpadlib==1.10.16
+lazr.restfulclient==0.14.4
+lazr.uri==1.0.6
+lightdm-gtk-greeter-settings==1.2.2
+Mako==1.1.3
+MarkupSafe==2.0.1
+matplotlib==3.6.1
+menulibre==2.2.2
+more-itertools==8.10.0
+mugshot==0.4.3
+netifaces==0.11.0
+numpy==1.23.4
+oauthlib==3.2.0
+olefile==0.46
+onboard==1.4.1
+packaging==21.3
+pexpect==4.8.0
+Pillow==9.0.1
+pluggy==1.0.0
+ply==3.11
+psutil==5.9.0
+ptyprocess==0.7.0
+pycairo==1.20.1
+pycparser==2.21
+pycups==2.0.1
+PyGObject==3.42.1
+PyJWT==2.3.0
+pymacaroons==0.13.0
+PyNaCl==1.5.0
+pyparsing==2.4.7
+pytest==7.2.1
+python-apt==2.3.0+ubuntu2.1
+python-dateutil==2.8.1
+python-debian===0.1.43ubuntu1
+PyYAML==5.4.1
+reportlab==3.6.8
+scikit-learn==1.1.2
+scipy==1.9.2
+SecretStorage==3.3.1
+sgt-launcher==0.2.7
+six==1.16.0
+sklearn==0.0
+ssh-import-id==5.11
+subprocess.run==0.0.8
+systemd-python==234
+threadpoolctl==3.1.0
+tomli==2.0.1
+ubuntu-advantage-tools==27.12
+ubuntu-drivers-common==0.0.0
+ufw==0.36.1
+unattended-upgrades==0.1
+urllib3==1.26.5
+wadllib==1.3.6
+xcffib==0.11.1
+xkit==0.0.0
+zipp==1.0.0

src/baseline.py

@@ -0,0 +1,78 @@
+from BaselineRemoval import BaselineRemoval 
+import numpy as np
+from scipy.spatial import ConvexHull
+from scipy import sparse
+
+# rubberband
+
+def baseline_zhang(y, polynomial_degree=2):
+    """
+	adaptive iteratively reweighted Penalized Least Squares (airPLS) - Zhi-Min Zhang et.al
+	https://pubs.rsc.org/is/content/articlelanding/2010/an/b922045c#!divAbstract
+	"""
+    baseObj = BaselineRemoval(y)
+    return baseObj.ZhangFit()
+
+def rubberband(x, y):
+    base = ConvexHull(list(zip(x, y))).vertices
+    base = np.roll(base, -base.argmax() - 1)
+    base1 = base[base.argmin():]
+    base2 = base[:base.argmin() + 1]
+    base1 = list(base1 if y[base1[1]] < y[base2[1]] else base2)
+    base1 =  [len(x) - 1] + base1 + [0]
+    new_y = y - np.interp(x, x[base1], y[base1])
+    return x, new_y
+
+def baseline_rubberband(y):
+    size = len(y)
+    im = [0]
+    for i in range(1, size - 1):
+        if y[i - 1] > y[i] < y[i + 1]:
+            im.append(i)
+    im.append(size - 1)
+    i_gm = np.argmin(y)  # index of global minimum
+    i_lm = i_gm  # index of local minimum
+    max_delta = np.max(y) - np.min(y)
+    indexies = []
+    while i_lm < size - 1:
+        min_slope = max_delta
+        i0 = im.index(i_lm)
+        i_lm = size - 1
+        for i in im[i0:]:
+            slope = (y[i] - y[im[i0]]) / float(i - im[i0])
+            if slope < min_slope:
+                min_slope = slope
+                i_lm = i
+        indexies.append(i_lm)
+    i_lm = i_gm
+    while i_lm > 1:
+        min_slope = max_delta
+        i0 = im.index(i_lm)
+        i_lm = 0
+        for i in im[i0:0:-1]:
+            slope = (y[i] - y[im[i0]]) / float(i - im[i0])
+            if slope < min_slope:
+                min_slope = slope
+                i_lm = i
+        indexies.append(i_lm)
+    indexies.append(i_gm)
+    indexies = sorted(list(set(indexies)))
+    z = np.interp(range(size), indexies, y[indexies])
+    return y - z
+
+def baseline_alss(y, lam=1e6, p=1e-3, niter=10):
+    """
+	an algorithm called "Asymmetric Least Squares Smoothing" by P. Eilers and H. Boelens
+	https://stackoverflow.com/questions/29156532/python-baseline-correction-library
+	"""
+    L = len(y)
+    D = sparse.diags([1, -2, 1], [0, -1, -2], shape=(L, L - 2))
+    D = lam * D.dot(D.transpose())
+    w = np.ones(L)
+    W = sparse.spdiags(w, 0, L, L)
+    for i in range(niter):
+        W.setdiag(w)
+        Z = W + D
+        z = spsolve(Z, w * y)
+        w = p * (y > z) + (1 - p) * (y < z)
+    return y - z

src/deconvolution.py

@@ -0,0 +1,175 @@
+import numpy as np
+import matplotlib.pyplot as plt
+from scipy.signal import savgol_filter
+from scipy.optimize import minimize
+from spectrum import Spectrum, show_spectra
+from enums import LossFunc
+from random import random
+
+
+def gauss(x, mu, sigma):
+    return 1./sigma/np.sqrt(2.*np.pi) * np.exp(-0.5 * ((x - mu) / sigma)**2)
+
+def lorentz(x, x0, gamma):
+    return 1./np.pi/gamma * (gamma**2 / (gamma**2 + (x - x0)**2))
+
+
+
+class Deconvolutor:
+
+    def __init__(self, spectrum, loss_func=LossFunc.RMSE, threshold=0.01):
+        self.spectrum = spectrum * 1
+        self.threshold = threshold
+        self.minima = None
+        self.loss_func = loss_func
+
+    @staticmethod
+    def __split_params(v):
+        assert len(v) % 3 == 0
+        assert len(v) != 0
+        third = len(v) // 3
+        ams, mus, sigmas = v[:third], v[third: 2 * third], v[2 * third:]
+        return ams, mus, sigmas
+
+    def loss(self, v):
+        # ams, mus, sigmas = Deconvolutor.__split_params(v)
+        # approx_data = np.zeros(len(data))
+        # for a, m, s in zip(ams, mus, sigmas):
+        #     approx_data += a * gauss(wavenums, m, s)
+        return self.loss_func(self.spectrum.data, self.__approximation(v)) # np.sum(np.square(data - approx_data))
+
+    def deconvolute(self):
+        #spc = Spectrum()
+        # plt.figure()
+        spc = self.spectrum * 1
+        # plt.plot(self.spectrum.wavenums, self.spectrum.data, 'b')
+        spc.data = spc.get_derivative(n=2)
+        # plt.plot(spc.wavenums, spc.data, 'p')
+        # plt.show()
+        indices = spc.get_extrema(minima=True)[0]
+        thr = spc.std * self.threshold
+        indices = filter(lambda x: abs(self.spectrum.data[x]) > thr, indices)
+        self.minima = list(indices)
+        # print('Indices: ', self.minima)
+        initial_params = self.__init_params()
+        res = minimize(self.loss, #lambda x: self.loss(x, data=self.spectrum.data),
+                       initial_params, method='Nelder-Mead', tol=1e-6).x
+        return zip(*Deconvolutor.__split_params(res))
+
+    def __approximation(self, params_vector):
+        approx_data = np.zeros(len(self.spectrum))
+        for a, m, s in zip(*Deconvolutor.__split_params(params_vector)):
+            approx_data += a * Deconvolutor.distribution()(self.spectrum.wavenums, m, s)
+        return approx_data
+
+    @classmethod
+    def distribution(cls, gauss_proba=1):
+        if random() < gauss_proba:
+            return gauss
+        else:
+            return lorentz
+
+    def __init_params(self):
+        amplitudes = np.ones(len(self.minima))
+        mus = np.array(self.spectrum.wavenums[self.minima])
+        sigmas = 1. / np.array(self.spectrum.data[self.minima]) / np.sqrt(2. * np.pi)
+        # amplitudes = self.spectrum.data[self.minima]
+        # mus = self.spectrum.wavenums[self.minima]
+        # sigmas = 1. / self.spectrum.data[self.minima] / np.sqrt(2. * np.pi)
+        return np.concatenate((amplitudes, mus, sigmas))
+
+# def split_params(v):
+#     assert len(v) % 3 == 0
+#     assert len(v) != 0
+#     third = len(v) // 3
+#     ams, mus, sigmas = v[:third], v[third: 2 * third], v[2 * third:]
+#     return ams, mus, sigmas
+
+#
+if __name__ == '__main__':
+    print('Hi')
+    wavenums = np.arange(4000., 600., -0.9367)
+    data = 1.5 * gauss(wavenums, 2567, 100) + gauss(wavenums, 1230, 200) * 2 + gauss(wavenums, 840, 29)
+    spc = Spectrum(wavenums, data)
+    dcv = Deconvolutor(spc)
+    approx = np.zeros(len(spc))
+    for a, m, s in dcv.deconvolute():
+        approx += a * gauss(spc.wavenums, m, s)
+    approx = Spectrum(spc.wavenums, approx, clss='app')
+    show_spectra([approx, spc])
+    # print(*[i for i in dcv.deconvolute()], sep='\n')
+
+
+# wavenums = np.arange(4000., 600., -0.9367)
+# data = 1.5 * gauss(wavenums, 2567, 100) + gauss(wavenums, 1230, 200) * 2 + gauss(wavenums, 840, 29)
+# spc= Spectrum(wavenums, data)
+#
+# def get_derivative(data, n=1, win_wight=13, order=5):
+# 	return savgol_filter(
+# 		data, win_wight, polyorder=order, deriv=n)
+#
+#
+# der2 = spc.get_derivative(n=2)
+#
+# minima = []
+# for i in range(1, len(der2) - 1):
+#     if der2[i-1] > der2[i] < der2[i+1]:
+#         minima.append(i)
+# print('Der2 minima', minima)
+# threshold = 0.01 * np.std(data)
+# minima_wn = []
+# minima_data = []
+# for i in minima:
+#     if data[i] >= threshold:
+#         minima_wn.append(wavenums[i])
+#         minima_data.append(data[i])
+# print('Threshold', threshold)
+# print('filtered wavenums: ', minima_wn)
+#
+# amplitudes = np.ones(len(minima_wn))
+# mus = np.array(minima_wn)
+# sigmas = 1./np.array(minima_data)/np.sqrt(2.*np.pi)
+# # minima_wn = np.array(minima_wn, dtype=int)
+# # amplitudes = data[minima_wn]
+# # mus = wavenums[minima_wn]
+# # sigmas = 1. / data[minima_wn] / np.sqrt(2. * np.pi)
+#
+# approx_data = np.zeros(len(data))
+# for a, m, s in zip(amplitudes, mus, sigmas):
+#     approx_data += a * gauss(wavenums, m, s)
+#
+# # plt.figure()
+# # plt.plot(wavenums, data, 'r', wavenums, approx_data,  'b')
+# # plt.show()
+#
+# def loss(v, data=data):
+#     third = len(v) // 3
+#     ams, mus, sigmas = v[:third], v[third: 2*third], v[2*third:]
+#     approx = np.zeros(len(data))
+#     for a, m, s in zip(ams, mus, sigmas):
+#         approx += a * gauss(wavenums, m, s)
+#     return np.sum(np.square(approx - data))
+#
+#
+# initial_params = np.array(list(amplitudes) + list(mus) + list(sigmas))
+# print("initial_params", initial_params)
+# res = minimize(lambda x: loss(x, data=data), initial_params, method='Nelder-Mead', tol=1e-6).x
+# third = len(res) // 3
+# ams, mus, sigmas = res[:third], res[third: 2 * third], res[2*third:]
+# approx_data = np.zeros(len(data))
+# for a, m, s in zip(ams, mus, sigmas):
+#     print(a, m, s)
+#     approx_data += a * gauss(wavenums, m, s)
+#
+# plt.plot(wavenums, data, 'b', wavenums, approx_data, 'r')
+# plt.show()
+
+
+
+
+# x0 = [1.3, 0.7, 0.8, 1.9, 1.2]
+
+# res = minimize(rosen, x0, method='Nelder-Mead', tol=1e-6)
+
+# res.x
+# array([ 1.,  1.,  1.,  1.,  1.])

src/enums.py

@@ -0,0 +1,22 @@
+from enum import Enum
+import numpy as np
+
+class NormMode(Enum):
+    VECTOR = 'vector'
+    AREA = 'area'
+    MINMAX = 'minmax'
+
+class BaseLineMode(Enum):
+    RB = 'rb'
+    ALSS = 'alls'
+    ZHANG = 'zhang'
+
+class Scale(Enum):
+    WAVENUMBERS = 'Wavenumber cm-1'
+    WAVELENGTH_nm = 'Wavelength nm'
+    WAVELENGTH_um = 'Wavelength um'
+
+class LossFunc(Enum):
+    RMSE = lambda y, y0: np.sqrt(np.sum(np.square(y - y0)))
+    BINARY = lambda y, y0: np.sum(y == y0)
+    MAE = lambda y, y0: np.mean(np.abs(y - y0))