diff --git a/Figure_1.png b/Figure_1.png new file mode 100644 index 0000000..5033b40 Binary files /dev/null and b/Figure_1.png differ diff --git a/MCA/MainFrm.cpp b/MCA/MainFrm.cpp index 86999a2..15ceb26 100644 --- a/MCA/MainFrm.cpp +++ b/MCA/MainFrm.cpp @@ -368,8 +368,9 @@ void CMainFrame::OnOptOpen() { // 其他设置 m_bClearFlag = TRUE; - pView->m_sRangeMode = "Auto"; + OnRangeAuto(); pView->m_sSmoothType = "Origin"; + UpdateValue(); pDoc->UpdateAllViews(NULL); } @@ -519,8 +520,8 @@ void CMainFrame::OnAxisLog() { if (pView == nullptr) pView = (CControlView*)pDoc->GetView(RUNTIME_CLASS(CControlView)); pView->m_sAxisMode = (CString)"Log"; - pView->m_sRangeMode = "Auto"; pView->m_ComboAxis.SetCurSel(1); + OnRangeAuto(); m_bAutoFlag = FALSE; m_bLinearFlag = TRUE; m_bLogFlag = FALSE; @@ -538,6 +539,7 @@ void CMainFrame::OnRangeAuto() { pView->m_sRangeMode = (CString)"Auto"; m_bAutoFlag = FALSE; pView->m_nLC = max(1, RoundupPowerof2(pView->m_nMaxCount)); + pView->UpdateData(FALSE); pDoc->UpdateAllViews(NULL); } @@ -550,7 +552,7 @@ void CMainFrame::OnRangeD4() { if (pView == nullptr) pView = (CControlView*)pDoc->GetView(RUNTIME_CLASS(CControlView)); pView->m_sRangeMode = (CString)"Manual"; - pView->m_nLC /= 4; + pView->m_nLC = pView->m_nLC > 1 ? pView->m_nLC / 4 : 1; m_bAutoFlag = TRUE; pView->UpdateData(FALSE); pDoc->UpdateAllViews(NULL); @@ -561,7 +563,7 @@ void CMainFrame::OnRangeD2() { if (pView == nullptr) pView = (CControlView*)pDoc->GetView(RUNTIME_CLASS(CControlView)); pView->m_sRangeMode = (CString)"Manual"; - pView->m_nLC /= 2; + pView->m_nLC = pView->m_nLC > 1 ? pView->m_nLC / 2: 1; m_bAutoFlag = TRUE; pView->UpdateData(FALSE); pDoc->UpdateAllViews(NULL); diff --git a/MCA/Matrix.h b/MCA/Matrix.h index d1b7c22..872847e 100644 --- a/MCA/Matrix.h +++ b/MCA/Matrix.h @@ -1,10 +1,10 @@ #pragma once -#ifndef DESCSS_Matrix_h -#define DESCSS_Matrix_h +#ifndef Matrix_h +#define Matrix_h #include -#define swap2(x, y) ((x)=(x)^(y),(y)=(x)^(y),(x)=(x)^(y)) +#define swap2(x, y) ((x) = (x) ^ (y), (y) = (x) ^ (y), (x) = (x) ^ (y)) class matrix { public: @@ -71,15 +71,9 @@ void matrix::Input(double* p) { for (int j = 0; j < m; j++) this->element[i][j] = *(p++); } -void matrix::Output() { - for (int i = 0; i < n; i++) { - for (int j = 0; j < m; j++) std::cout << element[i][j] << " "; - std::cout << std::endl; - } -} - double& matrix::operator()(int i, int j) const { - if (i < 0 || i > n - 1 || j < 0 || j > m - 1) std::cout << "Matrix Index Out Of Bounds " << std::endl; + if (i < 0 || i > n - 1 || j < 0 || j > m - 1) + std::cout << "Matrix Index Out Of Bounds " << std::endl; return element[i][j]; } @@ -123,15 +117,22 @@ matrix matrix::operator/(const double& C) const { return w; } -BOOL LUPDescomposition(matrix mA, matrix& mL, matrix& mU, int* P) -{ +std::ostream& operator<<(std::ostream& cout, matrix& mA) { + for (int i = 0; i < mA.rows(); i++) { + for (int j = 0; j < mA.cols(); j++) cout << mA(i, j) << " "; + if (i < mA.rows() - 1) cout << std::endl; + } + return cout; +} + +bool LUPDescomposition(matrix mA, matrix& mL, matrix& mU, int* P) { int n = mA.rows(); int row = 0, tmp; double L, U, tmp2; - for (int i = 0; i < n; i++) P[i] = i; + for (int i = 0; i < n; i++) P[i] = i; for (int i = 0; i < n - 1; i++) { tmp2 = 0.; - for(int j = i; j < n; j++) + for (int j = i; j < n; j++) if (std::abs(mA(j, i)) > tmp2) { tmp2 = std::abs(mA(j, i)); row = j; @@ -144,23 +145,21 @@ BOOL LUPDescomposition(matrix mA, matrix& mL, matrix& mU, int* P) tmp = P[i], P[i] = P[row], P[row] = tmp; - for (int j = 0; j < n; j++) - tmp2 = mA(i, j), mA(i, j) = mA(row, j), mA(row, j) = tmp2; + for (int j = 0; j < n; j++) tmp2 = mA(i, j), mA(i, j) = mA(row, j), mA(row, j) = tmp2; U = mA(i, i), L = 0; for (int j = i + 1; j < n; j++) { L = mA(j, i) / U; mA(j, i) = L; - for (int k = i + 1; k < n; k++) - mA(j, k) = mA(j, k) - mA(i, k) * L; + for (int k = i + 1; k < n; k++) mA(j, k) = mA(j, k) - mA(i, k) * L; } } for (int i = 0; i < n; i++) { mL(i, i) = 1; - for (int j = 0; j < i; j++) mL(i, j) = mA(i, j); - for (int j = i; j < n; j++) mU(i, j) = mA(i, j); + for (int j = 0; j < i; j++) mL(i, j) = mA(i, j); + for (int j = i; j < n; j++) mU(i, j) = mA(i, j); } - return TRUE; + return true; } matrix LUPSolve(matrix& mL, matrix& mU, int* P, double* b) { @@ -170,12 +169,12 @@ matrix LUPSolve(matrix& mL, matrix& mU, int* P, double* b) { for (int i = 0; i < n; i++) { y[i] = b[P[i]]; - for (int j = 0; j < i; j++) y[i] -= mL(i, j) * y[j]; + for (int j = 0; j < i; j++) y[i] -= mL(i, j) * y[j]; } for (int i = n - 1; i >= 0; i--) { mX(i, 0) = y[i]; - for (int j = n - 1; j > i; j--) mX(i, 0) -= mU(i, j) * mX(j, 0); + for (int j = n - 1; j > i; j--) mX(i, 0) -= mU(i, j) * mX(j, 0); mX(i, 0) /= mU(i, i); } @@ -198,15 +197,15 @@ matrix matrix::I() { } for (int i = 0; i < n; i++) - for (int j = 0; j < m; j++) mMir.element[i][j] = element[i][j]; + for (int j = 0; j < m; j++) mMir.element[i][j] = element[i][j]; LUPDescomposition(mMir, mL, mU, P); for (int i = 0; i < n; i++) { - for (int j = 0; j < n; j++) b[j] = 0; + for (int j = 0; j < n; j++) b[j] = 0; b[i] = 1; mX = LUPSolve(mL, mU, P, b); - for (int j = 0; j < n; j++) mInv(j, i) = mX(j, 0); + for (int j = 0; j < n; j++) mInv(j, i) = mX(j, 0); } return mInv; } @@ -219,16 +218,17 @@ matrix matrix::T() { return w; } -double Gaussian(double x, double A, double mean, double sigma) { +double Gaussian_(double x, double A, double mean, double sigma) { return A * exp(-(x - mean) * (x - mean) / (2 * sigma * sigma)); } -matrix Jacobian(int n, double* X, double A, double mean, double sigma) { +matrix GaussianJacobian_(int n, double* X, double A, double mean, double sigma) { matrix J(n, 3); for (int i = 0; i < n; i++) { - J(i, 0) = -Gaussian(X[i], A, mean, sigma) / A; - J(i, 1) = -(X[i] - mean) * Gaussian(X[i], A, mean, sigma) / (sigma * sigma); - J(i, 2) = -(X[i] - mean) * (X[i] - mean) * Gaussian(X[i], A, mean, sigma) / (sigma * sigma * sigma); + J(i, 0) = -Gaussian_(X[i], A, mean, sigma) / A; + J(i, 1) = -(X[i] - mean) * Gaussian_(X[i], A, mean, sigma) / (sigma * sigma); + J(i, 2) = -(X[i] - mean) * (X[i] - mean) * Gaussian_(X[i], A, mean, sigma) / + (sigma * sigma * sigma); } return J; }