CreateMatchFilter in main

This commit is contained in:
kradchen
2023-04-28 13:22:53 +08:00
parent d9e4fc5a0d
commit 9037838b5e

View File

@@ -6,7 +6,6 @@
#include <complex> #include <complex>
#include "Matrix.h" #include "Matrix.h"
#include "Function.h" #include "Function.h"
#include "Function1D.h" #include "Function1D.h"
@@ -14,63 +13,92 @@
#include "Function3D.h" #include "Function3D.h"
#include "MatlabReader.h" #include "MatlabReader.h"
int main() { int main()
{
MatlabReader m; MatlabReader m;
Input i; Input i;
MatchedFilter o; MatchedFilter o;
m.read(&i,&o,"/home/krad/TestData/testCreateMatchedFilter.mat"); m.read(&i, &o, "/home/krad/TestData/testCreateMatchedFilter.mat");
bool measuredCEused = true,findDefects = i.mParams->mFindDefects,removeOutliers = i.mRemoveOutliersFromCEMeasured; bool measuredCEused = true, findDefects = i.mParams->mFindDefects, removeOutliers = i.mRemoveOutliersFromCEMeasured;
Aurora::Matrix mFTime = Aurora::Matrix::fromRawData(i.mCe,4000,2304); Aurora::Matrix mFTime = Aurora::Matrix::fromRawData(i.mCe, 4000, 2304);
if (removeOutliers){ if (removeOutliers)
{
auto normSTD = Aurora::std(Aurora::abs(mFTime)); auto normSTD = Aurora::std(Aurora::abs(mFTime));
Aurora::nantoval(normSTD,0.0); Aurora::nantoval(normSTD, 0.0);
auto sortSTD = Aurora::sort(normSTD); auto sortSTD = Aurora::sort(normSTD);
int t = (int)std::round(0.4*mFTime.getDimSize(1))-1; int t = (int)std::round(0.4 * mFTime.getDimSize(1)) - 1;
t = t<=0?1.0:t; t = t <= 0 ? 1.0 : t;
auto absFTime = abs(mFTime); auto absFTime = abs(mFTime);
auto maxAbsFTime = max(absFTime); auto maxAbsFTime = max(absFTime);
auto maxFlag = maxAbsFTime < (0.1 * max(maxAbsFTime)); auto maxFlag = maxAbsFTime < (0.1 * max(maxAbsFTime));
auto lessFlag = normSTD < sortSTD(0,t).toMatrix(); auto lessFlag = normSTD < sortSTD(0, t).toMatrix();
long maxCol,maxRow; long maxCol, maxRow;
max(normSTD,Aurora::Column,maxRow,maxCol); max(normSTD, Aurora::Column, maxRow, maxCol);
for (int j = 0; j < mFTime.getDimSize(1); ++j) { for (int j = 0; j < mFTime.getDimSize(1); ++j)
if((bool)(lessFlag.getData()[j])||(bool)(maxFlag.getData()[j])){ {
mFTime(Aurora::$,j) = mFTime(Aurora::$,maxCol); if ((bool)(lessFlag.getData()[j]) || (bool)(maxFlag.getData()[j]))
{
mFTime(Aurora::$, j) = mFTime(Aurora::$, maxCol);
} }
} }
} }
auto matchedFilter = fft(mFTime); auto matchedFilter = fft(mFTime);
auto sumDiff = Aurora::zeros(1, matchedFilter.getDimSize(1));
long minCol = 998, row;
auto absMatchedFilter = abs(matchedFilter); auto absMatchedFilter = abs(matchedFilter);
// auto highNoiseScore = mean(absMatchedFilter) * Aurora::std(absMatchedFilter); auto highNoiseScore = mean(absMatchedFilter) * Aurora::std(absMatchedFilter);
auto highNoiseScore = mean( abs(matchedFilter)) * Aurora::std(absMatchedFilter); // auto highNoiseScore = mean( abs(matchedFilter)) * Aurora::std(absMatchedFilter);
printf("highNoiseScore[998]: %f4, highNoiseScore[2053]:%f4\r\n",highNoiseScore.getData()[998],highNoiseScore.getData()[2053]); printf("highNoiseScore[998]: %f4, highNoiseScore[2053]:%f4\r\n", highNoiseScore.getData()[998], highNoiseScore.getData()[2053]);
auto medianNoise = median(highNoiseScore); auto medianNoise = median(highNoiseScore);
printf("median: %f4\r\n",medianNoise.getScalar()); printf("median: %f , should be: 1724817516.8468074\r\n", medianNoise.getScalar());
long minCol,row;
min(abs(highNoiseScore - median(highNoiseScore)), Aurora::Column, row, minCol); min(abs(highNoiseScore - median(highNoiseScore)), Aurora::Column, row, minCol);
// //
minCol = 998; // minCol = 998;
auto sumDiff = Aurora::zeros(1, matchedFilter.getDimSize(1));
auto maxMatchFilter = matchedFilter(Aurora::$, minCol).toMatrix(); auto maxMatchFilter = matchedFilter(Aurora::$, minCol).toMatrix();
for (int k = 0; k < matchedFilter.getDimSize(1); ++k) { for (int k = 0; k < matchedFilter.getDimSize(1); ++k)
sumDiff.getData()[k] = sum(abs(matchedFilter(Aurora::$, k).toMatrix()- maxMatchFilter)).getData()[0]; {
sumDiff.getData()[k] = sum(abs(matchedFilter(Aurora::$, k).toMatrix() - maxMatchFilter)).getData()[0];
} }
auto indexe = sumDiff > (mean(sumDiff)+2.596 * Aurora::std(sumDiff).getScalar()); // printf("meanSumDiff\r\n");
for (int l = 0; l < indexe.getDataSize(); ++l) { // auto meanSumDiff = mean(sumDiff);
if ((bool)indexe.getData()[l]){ // printf("meanSumDiff finish\r\n");
matchedFilter(Aurora::$,l) = matchedFilter(Aurora::$,minCol); {
double meansumDiff = mean(sumDiff).getScalar() ;
double stdSumDiff = Aurora::std(sumDiff).getScalar();
double sumDiffJ =meansumDiff + 2.596 * stdSumDiff;
// auto indexe = sumDiff > sumDiffJ;
printf("indexe finish\r\n");
for (int l = 0; l < sumDiff.getDataSize(); ++l)
{
if (sumDiff.getData()[l]> sumDiffJ)
{
matchedFilter(Aurora::$, l) = matchedFilter(Aurora::$, minCol);
}
} }
printf("\r\n");
} }
if (measuredCEused){
auto mFTime2 = ifft(-matchedFilter);
max(Aurora::abs(mFTime2)).printfShape(); if (measuredCEused)
mFTime2 = mFTime2/repmat(max(Aurora::abs(mFTime2)),mFTime2.getDimSize(0),1); {
mFTime2 = mFTime2/repmat(sum(Aurora::abs(mFTime2)),mFTime2.getDimSize(0),1); auto mFTime2 = real(ifft(-matchedFilter));
mFTime2 = mFTime2 / repmat(max(Aurora::abs(mFTime2)), mFTime2.getDimSize(0), 1);
mFTime2 = mFTime2 / repmat(sum(Aurora::abs(mFTime2)), mFTime2.getDimSize(0), 1);
matchedFilter = fft(mFTime2); matchedFilter = fft(mFTime2);
} }
printf("run end\r\n"); for (size_t i = 0; i < 10000; i += 500)
{
printf("index :%d,origin output:%f4,%f4, output:%f4,%f4\r\n",
i, o.mReal[i], o.mImag[i], matchedFilter.getData()[i * 2], matchedFilter.getData()[i * 2 + 1]);
}
return 0; return 0;
} }