From cd286c72f149a11e434b8e1b09ffc9e805ed2699 Mon Sep 17 00:00:00 2001 From: YiHui Liu Date: Mon, 11 Jul 2022 09:18:23 +0800 Subject: [PATCH] add: astropy;change: bind fit --- .gitignore | 3 ++ main.py | 19 +++++---- qdx/Bind.py | 103 ++++++++++++++++++++++++++++++++++++++++------ qdx/BindFilter.py | 4 +- qdx/Block.py | 1 - requirements.txt | Bin 140 -> 168 bytes 6 files changed, 107 insertions(+), 23 deletions(-) diff --git a/.gitignore b/.gitignore index 3c2fb15..0189a86 100644 --- a/.gitignore +++ b/.gitignore @@ -17,5 +17,8 @@ __pycache__ .vscode .vs +# dev +*.ipynb + *.code-workspace *.exe diff --git a/main.py b/main.py index 1a4e04c..b76f213 100644 --- a/main.py +++ b/main.py @@ -10,21 +10,24 @@ BF = BindFilter() path = 'result/bind' file_list = os.listdir(path) -reg = re.compile(r'(([0-9]{4})-([0-9])-([0.9])).txt') +reg = re.compile(r'(([0-9]{4})-([0-9])-([0-9])).txt') pbar = tqdm(desc="Gaussian Mixture Bind Filter", total=len(file_list)) - for file in file_list: name, E, n, m = reg.match(file).groups() file = os.path.join(path, file) - BF(file) - BF.draw('result/GMM/' + name + '.png') - np.savetxt('result/bind-GMM/' + name + '.txt', BF.fit_data, fmt='%d') + BF.filter(file) + # BF.draw('result/GMM/' + name + '.png') + # np.savetxt('result/bind-GMM/' + name + '.txt', BF.fit_data, fmt='%d') binds.append(Bind(int(E), int(n), int(m), BF.fit_data)) pbar.update(1) - - break - +pbar.close() + +pbar = tqdm(desc="Bind Linear Fit", total=len(binds)) +for bind in binds: + bind.fit() + bind.draw('result/L-FIT/' + bind.name + '.png') + pbar.update(1) pbar.close() diff --git a/qdx/Bind.py b/qdx/Bind.py index d502932..b64f5a6 100644 --- a/qdx/Bind.py +++ b/qdx/Bind.py @@ -1,28 +1,107 @@ -from cProfile import label import numpy as np import matplotlib.pyplot as plt +from astropy.modeling import models, fitting +from sklearn.mixture import GaussianMixture from sklearn.linear_model import LinearRegression +def get_hist(data, maxN=50): + step = 1 + edge = np.arange(data.min(), data.max() + 1, step) + count, _ = np.histogram(data, bins=edge) + while count.max() <= maxN: + step += 1 + edge = np.arange(data.min(), data.max() + 1, step) + count, _ = np.histogram(data, bins=edge) + return count, (edge[1:] + edge[:-1]) / 2 class Bind(object): - k, b = 0, 0 + flag = 0 + L, C = 40, 0 + k1, k2, b = 0, 0, 0 + K1, C1, K2, C2 = 0, 0, 0, 0 - def __init__(self, E, n, m, data): - self.E = E + def __init__(self, n, m): self.n, self.m = n, m - self.x, self.y = data[:, 0], data[:, 1] + + def add_data(self, E, data): + if self.flag == 0: + self.flag = 1 + self.E1 = E + self.x1, self.y1 = data[:, 0].reshape(-1, 1), data[:, 1].reshape(-1, 1) + else: + self.E2 = E + self.x2, self.y2 = data[:, 0].reshape(-1, 1), data[:, 1].reshape(-1, 1) def fit(self): - self.reg = LinearRegression() - self.reg.fit(self.x, self.y) - self.k = self.reg.coef_[0][0] - self.b = self.reg.intercept_[0] + self.reg1 = LinearRegression() + self.reg1.fit(self.x1, self.y1) + self.K1 = self.reg1.coef_[0][0] + self.C1 = self.reg1.intercept_[0] + + self.reg2 = LinearRegression() + self.reg2.fit(self.x2, self.y2) + self.K2 = self.reg2.coef_[0][0] + self.C2 = self.reg2.intercept_[0] + + def solve(self): + self.k1 = (self.E1 - self.E2) / (self.C1 - self.C2) + self.k2 = -(self.K1 + self.K2) / (self.k1 * 2) + self.b = (self.E1 - self.C1 * self.k1 + self.E2 - self.C2 * self.k1) / 2 + + def split(self, n=7): + self.cluster = [] + data = np.concatenate((self.x1, self.y1), axis=1) + + model = GaussianMixture(n_components=n) + model.fit(data) + + ny = model.predict(data) + for i in np.unique(ny): + idx = np.where(ny == i)[0] + self.cluster.append(data[idx]) + + self.cluster = np.array(self.cluster) + + def fit2(self): + self.fx = [] + self.fy = [] + fitter = fitting.LevMarLSQFitter() + + for data in self.cluster: + x1, x2 = data[:, 0], data[:, 1] + c1, e1 = get_hist(x1) + c2, e2 = get_hist(x2) + modelx1 = models.Gaussian1D(amplitude=c1.max(), mean=x1.mean(), stddev=x1.std()) + modelx2 = models.Gaussian1D(amplitude=c2.max(), mean=x2.mean(), stddev=x2.std()) + fitted_model1 = fitter(modelx1, e1, c1) + fitted_model2 = fitter(modelx2, e2, c2) + self.fx.append(fitted_model1.mean.value) + self.fy.append(fitted_model2.mean.value) + + def pX(self, x1=None, x2=None): + x1 = x1 if x1 else self.x1 + x2 = x2 if x2 else self.y1 + self.px = (self.k2 * x2 - self.k1 * x1) * self.L / self.E1 + self.C def draw(self, title): fig = plt.figure(figsize=(8, 8)) ax = fig.add_subplot(1, 1, 1) - ax.scatter(self.x, self.y, s=0.1, c='k', label='raw') - ax.plot(self.x, self.reg.predict(self.x), c='r', label='fit') + ax.scatter(self.x1, self.y1, s=0.1, c='black', label=r'$E_1$') + ax.scatter(self.x2, self.y2, s=0.1, c='dimgray', label=r'$E_2$') + ax.plot(self.x1, self.reg1.predict(self.x1), c='red', label=r'$x_2={:.4f}x_1+{:.4f},\ R^2={:.5f}$'.format(self.K1, self.C1, self.RSquare1)) + ax.plot(self.x2, self.reg2.predict(self.x2), c='orangered', label=r'$x_2={:.4f}x_1+{:.4f},\ R^2={:.5f}$'.format(self.K2, self.C2, self.RSquare2)) ax.legend() - fig.savefig(title) + fig.savefig(title, facecolor='w', transparent=False) plt.close() + + @property + def name(self): + return '{:d}-{:d}-{:d}-{:d}'.format(self.E1, self.E2, self.n, self.m) + + @property + def RSquare1(self): + return self.reg1.score(self.x1, self.y1) + + @property + def RSquare2(self): + return self.reg2.score(self.x2, self.y2) diff --git a/qdx/BindFilter.py b/qdx/BindFilter.py index c045261..16ee86f 100644 --- a/qdx/BindFilter.py +++ b/qdx/BindFilter.py @@ -2,12 +2,11 @@ import numpy as np from matplotlib import pyplot as plt from sklearn.mixture import GaussianMixture - class BindFilter(object): def __init__(self): pass - def fit(self, file): + def filter(self, file): self.clusters = [] self.fit_data = np.array([]) @@ -22,6 +21,7 @@ class BindFilter(object): self.fit_data = idx if len(idx) > len(self.fit_data) else self.fit_data self.clusters.append(data[idx]) + self.data = data self.fit_data = data[self.fit_data] def draw(self, title): diff --git a/qdx/Block.py b/qdx/Block.py index 53aee42..30c5f96 100644 --- a/qdx/Block.py +++ b/qdx/Block.py @@ -1,6 +1,5 @@ from .Bind import Bind - class Block(object): binds = [] diff --git a/requirements.txt b/requirements.txt index 2c74d8436865191074c0f44717ab592cdc5b811c..f3fe5eb50d8adaf61706355768507c17c4261e84 100644 GIT binary patch delta 38 scmeBST)`;!|6d|QF+&MM5ko#h0YfE&EfAVA=rI^F@G@{QOl0Z;0M3gDrT_o{ delta 10 RcmZ3%*u%*5|KCKpE&vvS1aANU