diff --git a/src/Aurora.cu b/src/Aurora.cu index eea20de..792a10c 100644 --- a/src/Aurora.cu +++ b/src/Aurora.cu @@ -66,23 +66,6 @@ __global__ void validKernel(const float* aData, const float* aValid, float* aOut } } -// __global__ void validSubKernel(const double* aValid, double* aOutput, unsigned int* aCount, int aValidSize) -// { -// int index = blockIdx.x * blockDim.x + threadIdx.x; -// if(index == 0) -// { -// for(int i=0;i tex; -cudaArray* array; - -__global__ void testKernel(float* aData,cudaTextureObject_t aTexObj, cudaSurfaceObject_t aSurface) -{ - float a = tex2D(tex,5.5,5.5); - float b = tex2D(aTexObj,5.5,5.5); - float2 c = tex2D(aSurface,1,1); - printf("%f\n",a); - printf("%f\n",b); - printf("%f\n",c.x); - printf("%f\n",c.y); -} - -__global__ void writeSurfaceKernel( cudaSurfaceObject_t aSurface) -{ - float2 value; - value.x = 100; - value.y = 99; - surf2Dwrite(value, aSurface, 1, 1 ); - -} - -void subTest(cudaTextureObject_t& aTexture) -{ - cudaResourceDesc resourceDesc; - cudaTextureDesc textureDesc; - memset(&resourceDesc, 0, sizeof(resourceDesc)); - resourceDesc.resType = cudaResourceTypeArray; - resourceDesc.res.array.array = array; // 指向设备端的 CUDA 数组 - -// 在 textureDesc 中设置纹理描述 - memset(&textureDesc, 0, sizeof(textureDesc)); - textureDesc.addressMode[0] = cudaAddressModeClamp; - textureDesc.addressMode[1] = cudaAddressModeClamp; - textureDesc.filterMode = cudaFilterModeLinear; - textureDesc.readMode = cudaReadModeElementType; - textureDesc.normalizedCoords = false; - //textureDesc.channelDesc = texChannelDescSpeedOfSoundField; - cudaCreateTextureObject(&aTexture, &resourceDesc, &textureDesc, nullptr); -} - -void Aurora::test(float* aData) -{ - tex.addressMode[0] = cudaAddressModeClamp; // Texturreferenz beschreiben - tex.addressMode[1] = cudaAddressModeClamp; - tex.filterMode = cudaFilterModeLinear; - tex.normalized = 0; - cudaChannelFormatDesc texChannelDescSpeedOfSoundField = cudaCreateChannelDesc(32, 0, 0, 0, cudaChannelFormatKindFloat); - cudaMallocArray(&array, &texChannelDescSpeedOfSoundField, 10, 9); - cudaMemcpyToArray(array, 0, 0, aData,10 * 9 *sizeof(float), cudaMemcpyHostToDevice); - cudaBindTextureToArray ( &tex, array, &texChannelDescSpeedOfSoundField ); - - cudaTextureObject_t textureObj; - subTest(textureObj); - - - - struct cudaResourceDesc resDesc; - memset(&resDesc, 0, sizeof(resDesc)); - resDesc.resType = cudaResourceTypeArray; - // Create the surface objects - resDesc.res.array.array = array; - cudaSurfaceObject_t inputSurfObj = 0; - cudaCreateSurfaceObject(&inputSurfObj, &resDesc); - - writeSurfaceKernel<<<1,1>>>(inputSurfObj); - cudaDeviceSynchronize(); - testKernel<<<1, 1>>>(aData,textureObj, inputSurfObj); - - - cudaDeviceSynchronize(); - cudaUnbindTexture(&tex); -} - void Aurora::sort(const Aurora::Matrix& aMatrix) { RECON_INFO("cuda start"); diff --git a/src/Aurora.h b/src/Aurora.h index e6de35b..995bcbc 100644 --- a/src/Aurora.h +++ b/src/Aurora.h @@ -6,18 +6,10 @@ #include "Matrix.h" namespace Aurora { - //__global__ void doubleToComplexKernel(const double* input, cufftDoubleComplex* output, int size); void doubleToComplex(const double* input, cufftDoubleComplex* output, int size); - - //__global__ void maxKernel(const float* aInput, const float* aOutput, int aSize); void max(const float* aInput, const float* aOutput, int aSize); - Aurora::CudaMatrix valid(const Aurora::CudaMatrix aData, const Aurora::CudaMatrix aValid); - void test(float* aData); - void sort(const Aurora::Matrix& aMatrix); - - //Aurora::CudaMatrix getTransmissionDataSubFunction(const Aurora::CudaMatrix& aFxMatrix, const Aurora::CudaMatrix& aFhMatrix); } diff --git a/src/common/convertfp16tofloat.cpp b/src/common/convertfp16tofloat.cpp deleted file mode 100644 index 39123c8..0000000 --- a/src/common/convertfp16tofloat.cpp +++ /dev/null @@ -1,107 +0,0 @@ -#include "convertfp16tofloat.h" - -#include "Function.h" -#include -#include -#include -namespace { - - // const ushort CONVERT_AND_VALUE = 15; - // // andblack - // const __m128i andBlock = _mm_set_epi16(15, 15, 15, 15, 15, 15, 15, 15); - // const __m128i andBlock2 = - // _mm_set_epi16(2047, 2047, 2047, 2047, 2047, 2047, 2047, 2047); - // const __m128i zeroBlock = _mm_set_epi16(0, 0, 0, 0, 0, 0, 0, 0); - // const __m128i oneBlock = _mm_set_epi16(1, 1, 1, 1, 1, 1, 1, 1); - // const __m128i twokBlock = - // _mm_set_epi16(2048, 2048, 2048, 2048, 2048, 2048, 2048, 2048); - // const uint CONVERT_ADD_VALUE = UINT32_MAX - 4095; - // void convert(short * ptr, double* des,bool single = false){ - // // 初始化值 - // auto value = _mm_set_epi16(ptr[0], ptr[1], ptr[2], ptr[3], single?ptr[0]:ptr[4], single?ptr[0]:ptr[5], - // single?ptr[0]:ptr[6], single?ptr[0]:ptr[7]); - // auto uvalue = _mm_set_epi16( - // (ushort)ptr[0], (ushort)ptr[1], (ushort)ptr[2], (ushort)ptr[3], - // (ushort)(single?ptr[0]:ptr[4]), (ushort)(single?ptr[0]:ptr[5]), - // (ushort)(single?ptr[0]:ptr[6]), (ushort)(single?ptr[0]:ptr[7])); - // // 位移 - // auto sign_bit = _mm_srli_epi16(value, 15); // 右移16位取符号位 - // auto exponent = _mm_srli_epi16(uvalue, 11); - // // and - // exponent = _mm_and_si128(exponent, andBlock); - // // and ，then convert to int 32 bits - // auto fraction3 = _mm256_cvtepi16_epi32(_mm_and_si128(uvalue, andBlock2)); - // auto hidden_bit_mask = - // (_mm_cmp_epi16_mask(sign_bit, oneBlock, _MM_CMPINT_EQ) & - // _mm_cmp_epi16_mask(exponent, zeroBlock, _MM_CMPINT_EQ)) | - // (_mm_cmp_epi16_mask(sign_bit, zeroBlock, _MM_CMPINT_EQ) & - // _mm_cmp_epi16_mask(exponent, zeroBlock, _MM_CMPINT_NE)); - // auto hidden_bit16 = _mm_maskz_set1_epi16(hidden_bit_mask, 2048); - // auto hidden_bit32 = _mm256_cvtepi16_epi32(hidden_bit16); - // auto outputBlock = _mm256_add_epi32(fraction3, hidden_bit32); - // auto sign_bit_add_value = _mm256_maskz_set1_epi32( - // _mm_cmp_epi16_mask(sign_bit, oneBlock, _MM_CMPINT_EQ), - // CONVERT_ADD_VALUE); - // outputBlock = _mm256_add_epi32(outputBlock, sign_bit_add_value); - // auto exponent_mask = - // _mm_cmp_epi16_mask(oneBlock, exponent, _MM_CMPINT_LT); - // exponent = _mm_sub_epi16(exponent, oneBlock); - // auto exponent32 = _mm256_cvtepi16_epi32(exponent); - // auto zeroBlock32 = _mm256_cvtepi16_epi32(zeroBlock); - // auto offsetCount = - // _mm256_mask_blend_epi32(exponent_mask, zeroBlock32, exponent32); - - - // outputBlock = _mm256_sllv_epi32(outputBlock, offsetCount); - - // des[3] = _mm256_extract_epi32(outputBlock, 4); - // des[2] = _mm256_extract_epi32(outputBlock, 5); - // des[1] = _mm256_extract_epi32(outputBlock, 6); - // des[0] = _mm256_extract_epi32(outputBlock, 7); - // if(single) return; - // des[7] = _mm256_extract_epi32(outputBlock, 0); - // des[6] = _mm256_extract_epi32(outputBlock, 1); - // des[5] = _mm256_extract_epi32(outputBlock, 2); - // des[4] = _mm256_extract_epi32(outputBlock, 3); - - // } -} - -Aurora::Matrix Recon::convertfp16tofloat(Aurora::Matrix aMatrix) { - - - // auto input = aMatrix.getData(); - // // uint16变换为float(32位)输出大小翻倍 - // auto output = Aurora::malloc(aMatrix.getDataSize() * 4); - // size_t rows = aMatrix.getDataSize() * sizeof(double) / sizeof(short); - // size_t total_count = aMatrix.getDataSize(); - - - // #pragma omp parallel for - // for (size_t i = 0; i < total_count; i += 8) { - // // 循环展开以避免过度的线程调用 - // if (i < total_count) { - // auto ptr = (short *)(input + i); - // double *des = output + i * 4; - // ::convert(ptr, des,i+1>total_count); - // } - // if (i+2 < total_count) { - // auto ptr = (short *)(input + i + 2); - // double *des = output + (i+2) * 4; - // ::convert(ptr, des,i+3>total_count); - // } - // if (i+4 < total_count) { - // auto ptr = (short *)(input + i + 4); - // double *des = output + (i+4) * 4; - // ::convert(ptr, des,i+5>total_count); - // } - // if (i+6 < total_count) { - // auto ptr = (short *)(input + i + 6); - // double *des = output + (i+6) * 4; - // ::convert(ptr, des,i+7>total_count); - // } - // } - // return Aurora::Matrix::New(output, aMatrix.getDimSize(0), - // aMatrix.getDimSize(1), aMatrix.getDimSize(2)); - -} \ No newline at end of file diff --git a/src/common/convertfp16tofloat.h b/src/common/convertfp16tofloat.h deleted file mode 100644 index f018169..0000000 --- a/src/common/convertfp16tofloat.h +++ /dev/null @@ -1,10 +0,0 @@ -#ifndef __CONVERTFP16TOFLOAT_H__ -#define __CONVERTFP16TOFLOAT_H__ -#include "Matrix.h" -namespace Recon { - - Aurora::Matrix convertfp16tofloat(Aurora::Matrix aMatrix); - -} - -#endif // __CONVERTFP16TOFLOAT_H__ \ No newline at end of file diff --git a/src/common/dataBlockCreation/getAScanBlockPreprocessed.cpp b/src/common/dataBlockCreation/getAScanBlockPreprocessed.cpp index f8bfb18..4f37937 100644 --- a/src/common/dataBlockCreation/getAScanBlockPreprocessed.cpp +++ b/src/common/dataBlockCreation/getAScanBlockPreprocessed.cpp @@ -1,6 +1,5 @@ #include "getAScanBlockPreprocessed.h" -#include "CudaMatrix.h" #include "Matrix.h" #include "blockingGeometryInfo.h" #include "removeDataFromArrays.h" @@ -11,36 +10,15 @@ #include "src/transmissionReconstruction/dataFilter/dataFilter.h" #include "src/reflectionReconstruction/dataFilter.h" -#include "Aurora.h" - using namespace Aurora; using namespace Recon; -#include -#include -void printTime() -{ - struct timeval tpend; - gettimeofday(&tpend,NULL); - int secofday = (tpend.tv_sec + 3600 * 8 ) % 86400; - int hours = secofday / 3600; - int minutes = (secofday - hours * 3600 ) / 60; - int seconds = secofday % 60; - int milliseconds = tpend.tv_usec/1000; - std::cout<< hours << ":" < 0) @@ -88,72 +62,6 @@ AscanBlockPreprocessed Recon::getAscanBlockPreprocessed(Parser* aParser, const A { result.ascanBlockPreprocessed = ascanBlock.ascanBlock; } -//printTime(); + return result; -} - - -AscanBlockPreprocessedCuda Recon::getAscanBlockPreprocessedCuda(Parser* aParser, const Aurora::Matrix& aMp, const Aurora::Matrix& aSl, const Aurora::Matrix& aSn, - const Aurora::Matrix& aRl, const Aurora::Matrix& aRn, GeometryInfo& aGeom, const MeasurementInfo& aMeasInfo, - bool aApplyFilter, bool aTransReco) -{ -//std::cout<<"strart"< 0) - { - result.ascanBlockPreprocessed = preprocessAscanBlockCuda(result.ascanBlockPreprocessed, aMeasInfo); - } - // else - // { - // result.ascanBlockPreprocessed = ascanBlock.ascanBlock; - // } -//printTime(); -//std::cout<<"end"< -namespace Recon { - bool notifyStart(const std::string& aReconID ); - bool notifyFinish(); - bool notifyProgress(int percent); -} -#endif // __NOTIFY_H__ \ No newline at end of file diff --git a/src/main.cxx b/src/main.cxx index 8ecd2d2..1971e62 100644 --- a/src/main.cxx +++ b/src/main.cxx @@ -1,19 +1,12 @@ #include "Parser.h" #include "Parser.h" #include "config/config.h" -#include "log/notify.h" -#include "log/notify.h" #include #include -#include "MatlabReader.h" #define EIGEN_USE_MKL_ALL #include "src/common/fileHelper.h" #include "startReconstructions.h" - -#include "transmissionReconstruction/detection/detection.h" -#include "transmissionReconstruction/detection/detection.cuh" -#include "startReconstructions.h" #include "log/log.h" /* 0 is data path. 1 is dataRef path. @@ -25,9 +18,9 @@ int main(int argc, char *argv[]) { int argNum = 5; std::vector args(argNum); - args[0] = "/home/AScans_Data/CAS0.1Phantom/20230823T165617/"; - args[1] = "/home/AScans_Data/CAS0.1Phantom/20230823T171851/"; - args[2] = "/home/krad/Storage/DICOM/00e04b741e9f_20240619T145752/"; + args[0] = "/home/sun/20230418T145123/"; + args[1] = "/home/sun/20230418T141000/"; + args[2] = "/home/sun/AscanData/"; args[3] = Recon::DEFAULT_CONFIG_PATH; argc = argc <= argNum? argc-1 : argNum; for (int i = 0; i < argc; i++) diff --git a/src/reflectionReconstruction/startReflectionReconstruction.cpp b/src/reflectionReconstruction/startReflectionReconstruction.cpp index 2afcf0b..97f84ec 100644 --- a/src/reflectionReconstruction/startReflectionReconstruction.cpp +++ b/src/reflectionReconstruction/startReflectionReconstruction.cpp @@ -9,7 +9,6 @@ #include "common/precalculateChannelList.h" #include "common/dataBlockCreation/getAScanBlockPreprocessed.h" #include "common/dataBlockCreation/removeDataFromArrays.h" -#include "log/notify.h" #include "reflectionReconstruction/preprocessData/determineOptimalPulse.h" #include "reflectionReconstruction/reconstructionSAFT/reconstructionSAFT.h" #include "src/reflectionReconstruction/preprocessData/preprocessAScanBlockForReflection.h" @@ -139,7 +138,6 @@ Aurora::Matrix Recon::startReflectionReconstruction( Parser* aParser, int aSAFT_ std::cout< @@ -33,7 +32,7 @@ using namespace Aurora; int Recon::startReconstructions( const std::string& aDataPath, const std::string& aDataRefPath, const std::string& aOutputPath) { - // MatlabWriter writer(aOutputPath); + MatlabWriter writer(aOutputPath + "URResult.mat"); Parser dataParser(aDataPath); Parser refParser(aDataRefPath); Recon::DICOMExporter exporter(dataParser.getPatientData(), dataParser.getMetaData()); @@ -96,7 +95,6 @@ int Recon::startReconstructions( const std::string& aDataPath, const std::string TempInfo tempRef; CEInfo ceRef; Matrix transformationMatricesRef; - //Recon::notifyProgress(1); if(transParams::runTransmissionReco) { expInfoRef = loadMeasurementInfos(&refParser); @@ -129,7 +127,6 @@ int Recon::startReconstructions( const std::string& aDataPath, const std::string transformationMatricesRef = Matrix(); motorPosAvailableRef = Matrix(); } - //Recon::notifyProgress(2); if(!ce.ce.isNull() && !ceRef.ce.isNull()) { Matrix isEqual = (ce.ce == ceRef.ce); @@ -156,12 +153,10 @@ int Recon::startReconstructions( const std::string& aDataPath, const std::string preComputes.matchedFilter = createMatchedFilter(ce.ceRef, ce.measuredCEUsed, reflectParams::findDefects, reconParams::removeOutliersFromCEMeasured, expInfo.Hardware); } } - //Recon::notifyProgress(3); if(expInfo.sampleRate != reflectParams::aScanReconstructionFrequency) { reflectParams::expectedAScanDataLength = ceil(expInfo.numberSamples * ((float)reflectParams::aScanReconstructionFrequency / expInfo.sampleRate)); } - //Recon::notifyProgress(4); TransmissionReconstructionResult transmissionResult; bool sosAvailable = false; bool attAvailable = false; @@ -220,26 +215,23 @@ int Recon::startReconstructions( const std::string& aDataPath, const std::string GeometryInfo geomRef = getGeometryInfo(motorPosAvailableRef, transformationMatricesRef, rlList, rnList, slList, snList); RECON_INFO("Start transmissionRecostruction."); - //Recon::notifyProgress(5); transmissionResult = startTransmissionReconstruction(mp_inter, mpRef_inter, slList_inter, snList_inter, rlList_inter, rnList_inter, temp, tempRef, geom, geomRef, expInfo, expInfoRef, preComputes, &dataParser, &refParser); attAvailable = true; sosAvailable = true; exporter.exportDICOM(transmissionResult.recoSOS, DICOMExporter::SOS); exporter.exportDICOM(transmissionResult.recoATT, DICOMExporter::ATT); - // writer.setMatrix(transmissionResult.recoSOS, "sos"); - // writer.setMatrix(transmissionResult.recoATT, "att"); + writer.setMatrix(transmissionResult.recoSOS, "sos"); + writer.setMatrix(transmissionResult.recoATT, "att"); } if(reflectParams::runReflectionReco) { Matrix recoSOS = transmissionResult.recoSOS; Matrix recoATT = transmissionResult.recoATT; - precalcImageParameters(geom); - //Recon::notifyProgress(21); + precalcImageParameters(geom); //检测可使用内存是否足够,输出警报用,todo //checkEnvAndMemory(reflectParams.imageInfos.IMAGE_XYZ); auto preProcessData = preprocessTransmissionReconstructionForReflection(recoSOS, recoATT, transmissionResult.ddmis, geom, temp); - //Recon::notifyProgress(22); Matrix mp_inter = intersect(motorPosAvailable, reflectParams::motorPos); Matrix slList_inter = intersect(slList, reflectParams::senderTasList); Matrix snList_inter = intersect(snList, reflectParams::senderElementList); @@ -252,14 +244,12 @@ int Recon::startReconstructions( const std::string& aDataPath, const std::string preComputes.offset = estimateOffset(expInfo, ce, preComputes.matchedFilter); reflectParams::gpuSelectionList = reconParams::gpuSelectionList; - //Recon::notifyProgress(25); RECON_INFO("Start reflectionRecostruction."); Matrix env = startReflectionReconstruction(&dataParser, preProcessData.saftMode, mp_inter, slList_inter, snList_inter, rlList_inter, rnList_inter, geom, preProcessData.transRecos, expInfo, preComputes); - // writer.setMatrix(env, "reflect"); + writer.setMatrix(env, "reflect"); exporter.exportDICOM(env, Recon::DICOMExporter::REFL); - //Recon::notifyProgress(99); } - // writer.write(); + writer.write(); return 0; } \ No newline at end of file diff --git a/src/transmissionReconstruction/dataPreperation.h b/src/transmissionReconstruction/dataPreperation.h index 707c95f..ba4b84e 100644 --- a/src/transmissionReconstruction/dataPreperation.h +++ b/src/transmissionReconstruction/dataPreperation.h @@ -3,8 +3,7 @@ #include "Matrix.h" namespace Recon { -Aurora::Matrix distanceBetweenTwoPoints(const Aurora::Matrix& aMPtsA, - const Aurora::Matrix& aMPtsB); +Aurora::Matrix distanceBetweenTwoPoints(const Aurora::Matrix& aMPtsA, const Aurora::Matrix& aMPtsB); Aurora::Matrix calculateWaterTemperature(Aurora::Matrix aMWaterTempS, Aurora::Matrix aMWaterTempR, diff --git a/src/transmissionReconstruction/detection/detection.cpp b/src/transmissionReconstruction/detection/detection.cpp index 22cff40..dedba67 100644 --- a/src/transmissionReconstruction/detection/detection.cpp +++ b/src/transmissionReconstruction/detection/detection.cpp @@ -13,7 +13,6 @@ #include "common/getMeasurementMetaData.h" #include "config/config.h" #include "calculateBankDetectAndHilbertTransformation.hpp" -#include "transmissionReconstruction/detection/detection.cuh" using namespace Aurora; namespace Recon { @@ -415,39 +414,34 @@ namespace Recon { return result; } - DetectResult transmissionDetection(const Aurora::CudaMatrix &AscanBlock, - const Aurora::CudaMatrix &AscanRefBlock, - const Aurora::CudaMatrix &distBlock, - const Aurora::CudaMatrix &distRefBlock, + DetectResult transmissionDetection(const Aurora::Matrix &AscanBlock, + const Aurora::Matrix &AscanRefBlock, + const Aurora::Matrix &distBlock, + const Aurora::Matrix &distRefBlock, const Aurora::Matrix &sosWaterBlock, const Aurora::Matrix &sosWaterRefBlock, float expectedSOSWater) { - auto _sosWaterBlock = temperatureToSoundSpeed(sosWaterBlock, "marczak").toDeviceMatrix(); - auto _sosWaterRefBlock = temperatureToSoundSpeed(sosWaterRefBlock, "marczak").toDeviceMatrix(); + auto _sosWaterBlock = temperatureToSoundSpeed(sosWaterBlock, "marczak"); + auto _sosWaterRefBlock = temperatureToSoundSpeed(sosWaterRefBlock, "marczak"); switch (Recon::transParams::version) { - // case 1: { - // return detectTofAndAttMex( - // AscanBlock, AscanRefBlock, distBlock, distRefBlock, - // _sosWaterBlock, _sosWaterRefBlock, Recon::transParams::resampleFactor, Recon::transParams::nThreads, - // expectedSOSWater, Recon::transParams::useTimeWindowing, - // Recon::transParams::aScanReconstructionFrequency, Recon::transParams::detectionWindowATT, - // Recon::transParams::offsetElectronic, Recon::transParams::detectionWindowSOS, Recon::transParams::minSpeedOfSound, - // Recon::transParams::maxSpeedOfSound, Recon::transParams::gaussWindow); - // } - // case 2: - default: - auto r = detectTofAndAtt( + case 1: { + return detectTofAndAttMex( + AscanBlock, AscanRefBlock, distBlock, distRefBlock, + _sosWaterBlock, _sosWaterRefBlock, Recon::transParams::resampleFactor, Recon::transParams::nThreads, + expectedSOSWater, Recon::transParams::useTimeWindowing, + Recon::transParams::aScanReconstructionFrequency, Recon::transParams::detectionWindowATT, + Recon::transParams::offsetElectronic, Recon::transParams::detectionWindowSOS, Recon::transParams::minSpeedOfSound, + Recon::transParams::maxSpeedOfSound, Recon::transParams::gaussWindow); + } + case 2: + default: + return detectTofAndAtt( AscanBlock, AscanRefBlock, distBlock, distRefBlock, _sosWaterBlock, _sosWaterRefBlock, Recon::transParams::resampleFactor, Recon::transParams::nThreads, expectedSOSWater, Recon::transParams::useTimeWindowing, Recon::transParams::aScanReconstructionFrequency, Recon::transParams::detectionWindowATT, Recon::transParams::offsetElectronic, Recon::transParams::detectionWindowSOS, Recon::transParams::minSpeedOfSound, Recon::transParams::maxSpeedOfSound, Recon::transParams::gaussWindow); - DetectResult ret; - ret.att = r.att.toHostMatrix(); - ret.sosValue = r.sosValue.toHostMatrix(); - ret.tof = r.tof.toHostMatrix(); - return ret; } } } diff --git a/src/transmissionReconstruction/detection/detection.h b/src/transmissionReconstruction/detection/detection.h index 76a5292..b30a360 100644 --- a/src/transmissionReconstruction/detection/detection.h +++ b/src/transmissionReconstruction/detection/detection.h @@ -75,8 +75,8 @@ detectTofAndAttMex( DetectResult transmissionDetection( - const Aurora::CudaMatrix &AscanBlock, const Aurora::CudaMatrix &AscanRefBlock, - const Aurora::CudaMatrix &distBlock, const Aurora::CudaMatrix &distRefBlock, + const Aurora::Matrix &AscanBlock, const Aurora::Matrix &AscanRefBlock, + const Aurora::Matrix &distBlock, const Aurora::Matrix &distRefBlock, const Aurora::Matrix &sosWaterBlock, const Aurora::Matrix &sosWaterRefBlock, float expectedSOSWater); } // namespace Recon diff --git a/src/transmissionReconstruction/detection/getTransmissionData.cpp b/src/transmissionReconstruction/detection/getTransmissionData.cpp index b78718e..2d206e6 100644 --- a/src/transmissionReconstruction/detection/getTransmissionData.cpp +++ b/src/transmissionReconstruction/detection/getTransmissionData.cpp @@ -1,12 +1,8 @@ #include "getTransmissionData.h" -#include "getTransmissionData.cuh" -#include "AuroraDefs.h" -#include "CudaMatrix.h" #include "Function.h" #include "Function1D.h" #include "Function2D.h" #include "common/dataBlockCreation/removeDataFromArrays.h" -#include "log/notify.h" #include "src//config/config.h" #include "src/common/getGeometryInfo.h" #include "src/common/temperatureCalculation/extractTasTemperature.h" @@ -31,11 +27,6 @@ #include #include -#include -#include "Function1D.cuh" -#include "Function2D.cuh" -#include - using namespace Recon; using namespace Aurora; @@ -50,8 +41,8 @@ namespace Matrix waterTempBlock; MetaInfos metaInfos; - CudaMatrix ascanBlock; - CudaMatrix ascanBlockRef; + Matrix ascanBlock; + Matrix ascanBlockRef; Matrix dists; Matrix distRefBlock; Matrix waterTempRefBlock; @@ -66,57 +57,26 @@ namespace int BUFFER_COUNT = 0; int BUFFER_SIZE = 4;//<=8 -void printTime() -{ - struct timeval tpend; - gettimeofday(&tpend,NULL); - int secofday = (tpend.tv_sec + 3600 * 8 ) % 86400; - int hours = secofday / 3600; - int minutes = (secofday - hours * 3600 ) / 60; - int seconds = secofday % 60; - int milliseconds = tpend.tv_usec/1000; - std::cout<< hours << ":" < lockBufferCount(CREATE_BUFFER_MUTEX); BLOCK_OF_TRANSIMISSIONDARA_BUFFER.erase(std::to_string(index)); @@ -392,7 +362,6 @@ TransmissionData Recon::getTransmissionData(const Aurora::Matrix& aMotorPos, con std::cout<<"Remove: "< #include #include #include -#include -#include using namespace Aurora; -using solveParameterIteratorFunctionType = std::vector> (*)(Aurora::Sparse M, Aurora::Matrix &b, - const Aurora::Matrix &dims, bool oneIter, bool nonNeg, int aDevice); -using slownessToSOSFunctionType = Matrix (*)(Aurora::Matrix & aVF1, float aSOS_IN_WATER); namespace Recon { Aurora::Matrix calculateMinimalMaximalTransducerPositions( const Aurora::Matrix &aMSenderList, const Aurora::Matrix &aMReceiverList) { @@ -195,9 +188,7 @@ namespace Recon { BuildMatrixResult buildMatrixR; for(int iter=1; iter<=numIter; ++iter) { - auto resDevice = res.toDeviceMatrix(); - //1200ms - buildMatrixR = buildMatrix(senderList.toDeviceMatrix(), receiverList.toDeviceMatrix(), resDevice, dims.toDeviceMatrix(), bentRecon && (iter!=1), potentialMap.toDeviceMatrix()); + buildMatrixR = buildMatrix(senderList, receiverList, res, dims, bentRecon && (iter!=1), potentialMap); if(!data.isNull() && bentRecon && iter != numIter) { //与默认配置bentRecon不符，暂不实现todo @@ -230,23 +221,25 @@ namespace Recon { { allHitMaps.push_back(buildMatrixR.hitmap); } - - if(!data.isNull()) - { - //1500ms - Matrix sosValue = solveParameterIterator(buildMatrixR.M, b, dims, false, transParams::nonNeg)[0][0]; - //1ms - result.outSOS = slownessToSOS(sosValue, SOS_IN_WATER) ; + #pragma omp parallel for num_threads(2) + for (int i =0; i<2; i++){ + if (i ==0){ + if(!data.isNull()) + { + Matrix sosValue = solveParameterIterator(buildMatrixR.M, b, dims, false, transParams::nonNeg)[0][0]; + result.outSOS = slownessToSOS(sosValue, SOS_IN_WATER) ; + } + } + else{ + if(!dataAtt.isNull()) + { + Matrix attValue = solveParameterIterator(buildMatrixR.M, bAtt, dims, false, transParams::nonNeg,1)[0][0]; + result.outATT = attValue/100 ; + } + } } + } - if(!dataAtt.isNull()) - { - //1500ms - Matrix attValue = solveParameterIterator(buildMatrixR.M, bAtt, dims, false, transParams::nonNeg)[0][0]; - //1ms - result.outATT = attValue/100 ; - } - } return result; } diff --git a/src/transmissionReconstruction/startTransmissionReconstruction.cpp b/src/transmissionReconstruction/startTransmissionReconstruction.cpp index cef6c53..b49ed73 100644 --- a/src/transmissionReconstruction/startTransmissionReconstruction.cpp +++ b/src/transmissionReconstruction/startTransmissionReconstruction.cpp @@ -1,10 +1,8 @@ #include "startTransmissionReconstruction.h" #include "./detection/getTransmissionData.h" #include "Matrix.h" -#include "CudaMatrix.h" #include "log/log.h" #include "common/dataBlockCreation/removeDataFromArrays.h" -#include "log/notify.h" #include "src/transmissionReconstruction/dataFilter/dataFilter.h" #include "src/transmissionReconstruction/dataPreperation.h" #include "src/common/getMeasurementMetaData.h" @@ -29,23 +27,19 @@ TransmissionReconstructionResult Recon::startTransmissionReconstruction(const Au auto transmissionData = getTransmissionData(aMotorPos, aMotoPosRef, aSlList, aSnList, aRlList, aRnList, aTemp, aTempRef, aGeom, aGeomRef, aExpInfo, aExpInfoRef, aPreComputes, aParser, aParserRef); Matrix dists = Recon::distanceBetweenTwoPoints(transmissionData.senderList, transmissionData.receiverList); - Matrix sosRef = Recon::temperatureToSoundSpeed(transmissionData.waterTempList, "marczak"); - //Recon::notifyProgress(17); + Matrix sosRef = Recon::temperatureToSoundSpeed(transmissionData.waterTempList, "marczak"); Matrix valid = Recon::checkTofDetections(transmissionData.tofDataTotal, dists, sosRef, - Recon::transParams::minSpeedOfSound,Recon::transParams::maxSpeedOfSound).valid; - //Recon::notifyProgress(18); + Recon::transParams::minSpeedOfSound,Recon::transParams::maxSpeedOfSound).valid; if(transParams::qualityCheck) { qualityReview(sum(valid,Aurora::All)[0], transmissionData.dataInfo.numPossibleScans); - } - //Recon::notifyProgress(19); + } DiscretizePositionValues positionValues = Recon::discretizePositions(transmissionData.senderList, transmissionData.receiverList, Recon::transParams::numPixelXY); Matrix tofData = removeDataFromArrays(transmissionData.tofDataTotal, valid); Matrix attData = removeDataFromArrays(transmissionData.attDataTotal, valid); Matrix senderList = removeDataFromArrays(positionValues.senderCoordList, valid); Matrix reveiverList = removeDataFromArrays(positionValues.receiverCoordList, valid); - //Recon::notifyProgress(20); RECON_INFO("Start reconstructArt."); auto transmissionReon = reconstructArt(tofData, attData, positionValues.dims, senderList, reveiverList, positionValues.res, aTemp.expectedSOSWater[0]); diff --git a/test/Common_Test.cpp b/test/Common_Test.cpp index 0457dd2..a75dc65 100644 --- a/test/Common_Test.cpp +++ b/test/Common_Test.cpp @@ -4,7 +4,6 @@ #include "Matrix.h" #include "common/ceMatchedFilterHandling.h" #include "common/common.h" -#include "common/convertfp16tofloat.h" #include "common/getGeometryInfo.h" #include "common/dataBlockCreation/getAscanBlock.h" #include "common/dataBlockCreation/blockingGeometryInfo.h" @@ -52,21 +51,6 @@ TEST_F(Common_Test, adaptFrequency) { } -TEST_F(Common_Test, convertfp16tofloat) { - MatlabReader m("/home/krad/TestData/convertReal.mat"); - - size_t count = 0; - auto input = m.readint16("input",count); - auto ma = Aurora::Matrix::copyFromRawData((float*)input.get(),count/4); - auto resultM = Recon::convertfp16tofloat(ma); - auto result = resultM.getData(); - auto output = m.read("output"); - for (size_t i = 0; i