diff --git a/Modules/OpenCL/mitkOclDataSet.cpp b/Modules/OpenCL/mitkOclDataSet.cpp index 514b2d100d..eaebdf5d42 100644 --- a/Modules/OpenCL/mitkOclDataSet.cpp +++ b/Modules/OpenCL/mitkOclDataSet.cpp @@ -1,160 +1,160 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #include "mitkOclDataSet.h" #include "mitkCommon.h" #include "mitkLogMacros.h" #include "mitkOclUtils.h" #include mitk::OclDataSet::OclDataSet() : m_gpuBuffer(nullptr), m_context(nullptr), m_bufferSize(0), m_gpuModified(false), m_cpuModified(false), m_Data(nullptr), m_BpE(1) { } mitk::OclDataSet::~OclDataSet() { MITK_INFO << "OclDataSet Destructor"; //release GMEM Image buffer if (m_gpuBuffer) clReleaseMemObject(m_gpuBuffer); } cl_mem mitk::OclDataSet::CreateGPUBuffer() { MITK_INFO << "InitializeGPUBuffer call with: BPE=" << m_BpE; us::ServiceReference ref = GetModuleContext()->GetServiceReference(); OclResourceService* resources = GetModuleContext()->GetService(ref); m_context = resources->GetContext(); int clErr; if (m_gpuBuffer) clReleaseMemObject(m_gpuBuffer); m_gpuBuffer = clCreateBuffer(m_context, CL_MEM_READ_WRITE, m_bufferSize * m_BpE, nullptr, &clErr); MITK_INFO << "Created GPU Buffer Object of size: " << m_BpE* m_bufferSize << " Bytes"; CHECK_OCL_ERR(clErr); return m_gpuBuffer; } bool mitk::OclDataSet::IsModified(int _type) { if (_type) return m_cpuModified; else return m_gpuModified; } void mitk::OclDataSet::Modified(int _type) { // defines... GPU: 0, CPU: 1 m_cpuModified = _type; m_gpuModified = !_type; } int mitk::OclDataSet::TransferDataToGPU(cl_command_queue gpuComQueue) { cl_int clErr = 0; // check whether an image present if (m_Data == nullptr){ MITK_ERROR("ocl.DataSet") << "(mitk) No data present!\n"; return -1; } // there is a need for copy only if RAM-Data newer then GMEM data if (m_cpuModified) { //check the buffer if(m_gpuBuffer == nullptr) { CreateGPUBuffer(); } if (m_gpuBuffer != nullptr) { clErr = clEnqueueWriteBuffer(gpuComQueue, m_gpuBuffer, CL_TRUE, 0, m_bufferSize * m_BpE, m_Data, 0, NULL, NULL); MITK_INFO << "Wrote Data to GPU Buffer Object."; } else { MITK_ERROR << "No GPU buffer present!"; } CHECK_OCL_ERR(clErr); m_gpuModified = true; } return clErr; } cl_mem mitk::OclDataSet::GetGPUBuffer() { // clGetMemObjectInfo() cl_mem_object_type memInfo; cl_int clErr = 0; // query image object info only if already initialized if( this->m_gpuBuffer ) { clErr = clGetMemObjectInfo(this->m_gpuBuffer, CL_MEM_TYPE, sizeof(cl_mem_object_type), &memInfo, nullptr ); CHECK_OCL_ERR(clErr); } MITK_INFO << "Querying info for object, recieving: " << memInfo; return m_gpuBuffer; } void* mitk::OclDataSet::TransferDataToCPU(cl_command_queue gpuComQueue) { cl_int clErr = 0; // if image created on GPU, needs to create mitk::Image if( m_gpuBuffer == nullptr ){ MITK_ERROR("ocl.DataSet") << "(mitk) No buffer present!\n"; return nullptr; } // check buffersize char* data = new char[m_bufferSize * m_BpE]; // debug info oclPrintMemObjectInfo( m_gpuBuffer ); - clErr = clEnqueueReadBuffer( gpuComQueue, m_gpuBuffer, CL_FALSE, 0, m_bufferSize * m_BpE, data ,0, nullptr, nullptr); + clErr = clEnqueueReadBuffer( gpuComQueue, m_gpuBuffer, CL_TRUE, 0, m_bufferSize * m_BpE, data ,0, nullptr, nullptr); CHECK_OCL_ERR(clErr); clFlush( gpuComQueue ); // the cpu data is same as gpu this->m_gpuModified = false; return (void*) data; } void mitk::OclDataSet::SetBufferSize(unsigned int size) { m_bufferSize = size; } void mitk::OclDataSet::SetBpE(unsigned short BpE) { m_BpE = BpE; } \ No newline at end of file diff --git a/Modules/PhotoacousticsAlgorithms/Resources/DAS.cl b/Modules/PhotoacousticsAlgorithms/Resources/DAS.cl index 2e8812e0dd..813e819292 100644 --- a/Modules/PhotoacousticsAlgorithms/Resources/DAS.cl +++ b/Modules/PhotoacousticsAlgorithms/Resources/DAS.cl @@ -1,65 +1,64 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ __kernel void ckDAS( __global float* dSource, // input image __global float* dDest, // output buffer __global unsigned short* usedLines, - __global unsigned int* memoryLocations, - __global unsigned short* AddSamples, + __global unsigned short* delays, __constant float* apodArray, unsigned short apodArraySize, unsigned int inputL, unsigned int inputS, unsigned int Slices, unsigned int outputL, unsigned int outputS // parameters ) { // get thread identifier unsigned int globalPosX = get_global_id(0); unsigned int globalPosY = get_global_id(1); unsigned int globalPosZ = get_global_id(2); // terminate non-valid threads if ( globalPosX < outputL && globalPosY < outputS && globalPosZ < Slices ) { + float l_i = (float)globalPosX / (float)outputL * (float)inputL; + unsigned short curUsedLines = usedLines[globalPosY * 3 * outputL + 3 * globalPosX]; unsigned short minLine = usedLines[globalPosY * 3 * outputL + 3 * globalPosX + 1]; unsigned short maxLine = usedLines[globalPosY * 3 *outputL + 3 * globalPosX + 2]; float apod_mult = (float)apodArraySize / (float)curUsedLines; - unsigned short AddSample = 0; + unsigned short Delay = 0; float output = 0; float mult = 0; - unsigned int MemoryStartAccessPoint = memoryLocations[globalPosY * outputL + globalPosX]; - for (short l_s = minLine; l_s < maxLine; ++l_s) { - AddSample = AddSamples[MemoryStartAccessPoint + l_s - minLine]; - if (AddSample < inputS && AddSample >= 0) { - output += apodArray[(int)((l_s - minLine)*apod_mult)] * dSource[(int)(globalPosZ * inputL * inputS + AddSample * inputL + l_s)]; + Delay = delays[globalPosY * inputL + (int)fabs(l_s - l_i)]; + if (Delay < inputS && Delay >= 0) { + output += apodArray[(int)((l_s - minLine)*apod_mult)] * dSource[(int)(globalPosZ * inputL * inputS + Delay * inputL + l_s)]; } else --curUsedLines; } dDest[ globalPosZ * outputL * outputS + globalPosY * outputL + globalPosX ] = output / (float)curUsedLines; } } \ No newline at end of file diff --git a/Modules/PhotoacousticsAlgorithms/Resources/DMAS.cl b/Modules/PhotoacousticsAlgorithms/Resources/DMAS.cl index f2ddcd5088..dfb1b5d334 100644 --- a/Modules/PhotoacousticsAlgorithms/Resources/DMAS.cl +++ b/Modules/PhotoacousticsAlgorithms/Resources/DMAS.cl @@ -1,77 +1,76 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ __kernel void ckDMAS( __global float* dSource, // input image __global float* dDest, // output buffer __global unsigned short* usedLines, - __global unsigned int* memoryLocations, __global unsigned short* AddSamples, __constant float* apodArray, unsigned short apodArraySize, unsigned int inputL, unsigned int inputS, unsigned int Slices, unsigned int outputL, unsigned int outputS // parameters ) { // get thread identifier unsigned int globalPosX = get_global_id(0); unsigned int globalPosY = get_global_id(1); unsigned int globalPosZ = get_global_id(2); // terminate non-valid threads if ( globalPosX < outputL && globalPosY < outputS && globalPosZ < Slices ) { + float l_i = (float)globalPosX / (float)outputL * (float)inputL; + unsigned short curUsedLines = usedLines[globalPosY * 3 * outputL + 3 * globalPosX]; unsigned short minLine = usedLines[globalPosY * 3 * outputL + 3 * globalPosX + 1]; unsigned short maxLine = usedLines[globalPosY * 3 *outputL + 3 * globalPosX + 2]; - + float apod_mult = (float)apodArraySize / (float)curUsedLines; - - unsigned short AddSample1 = 0; - unsigned short AddSample2 = 0; - + + unsigned short Delay1 = 0; + unsigned short Delay2 = 0; + float output = 0; float mult = 0; - - unsigned int MemoryStartAccessPoint = memoryLocations[globalPosY * outputL + globalPosX]; for (short l_s1 = minLine; l_s1 < maxLine; ++l_s1) { - AddSample1 = AddSamples[MemoryStartAccessPoint + l_s1 - minLine]; - if (AddSample1 < inputS && AddSample1 >= 0) { + Delay1 = AddSamples[globalPosY * inputL + (int)fabs(l_s1 - l_i)]; + if (Delay1 < inputS && Delay1 >= 0) { for (short l_s2 = l_s1 + 1; l_s2 < maxLine; ++l_s2) { - AddSample2 = AddSamples[MemoryStartAccessPoint + l_s2 - minLine]; - if (AddSample1 < inputS && AddSample1 >= 0) { + Delay2 = AddSamples[globalPosY * inputL + (int)fabs(l_s2 - l_i)]; + if (Delay1 < inputS && Delay1 >= 0) { mult = apodArray[(int)((l_s2 - minLine)*apod_mult)] * - dSource[(int)(globalPosZ * inputL * inputS + AddSample2 * inputL + l_s2)] + dSource[(int)(globalPosZ * inputL * inputS + Delay2 * inputL + l_s2)] * apodArray[(int)((l_s1 - minLine)*apod_mult)] * - dSource[(int)(globalPosZ * inputL * inputS + AddSample1 * inputL + l_s1)]; + dSource[(int)(globalPosZ * inputL * inputS + Delay1 * inputL + l_s1)]; output += sqrt(mult * ((float)(mult>0)-(float)(mult<0))) * ((mult > 0) - (mult < 0)); } } } else --curUsedLines; } - + dDest[ globalPosZ * outputL * outputS + globalPosY * outputL + globalPosX ] = output / (pow((float)curUsedLines, 2.0f) - (curUsedLines - 1)); } } \ No newline at end of file diff --git a/Modules/PhotoacousticsAlgorithms/Resources/DelayCalculation.cl b/Modules/PhotoacousticsAlgorithms/Resources/DelayCalculation.cl index 4f96c2f836..83649b305e 100644 --- a/Modules/PhotoacousticsAlgorithms/Resources/DelayCalculation.cl +++ b/Modules/PhotoacousticsAlgorithms/Resources/DelayCalculation.cl @@ -1,74 +1,72 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ __kernel void ckDelayCalculationQuad( __global unsigned short *gDest, - __global unsigned short *usedLines, - __global unsigned int *memoryLocations, + __global unsigned short *usedLines, unsigned int inputL, unsigned int inputS, unsigned int outputL, unsigned int outputS, char isPAImage, float delayMultiplicatorRaw // parameters ) - { - uint globalPosX = get_global_id(0); - uint globalPosY = get_global_id(1); - uint globalPosZ = get_global_id(2); - - if(globalPosZ < usedLines[globalPosY * 3 * outputL + 3 * globalPosX]) - { - float l_i = (float)globalPosX / outputL * inputL; - float s_i = (float)globalPosY / outputS * inputS / 2; - - float l_s = usedLines[globalPosY * 3 * outputL + 3 * globalPosX + 1] + globalPosZ; +{ + uint globalPosX = get_global_id(0); + uint globalPosY = get_global_id(1); + uint globalPosZ = get_global_id(2); + + if (globalPosX < inputL && globalPosY < outputS) + { + float l_i = 0; // we calculate the delays relative to line zero + float s_i = (float)globalPosY / (float)outputS * (float)inputS / 2; - float delayMultiplicator = delayMultiplicatorRaw / s_i; - float AddSample = delayMultiplicator * pow((l_s - l_i), 2) + s_i + (1-isPAImage)*s_i; - gDest[memoryLocations[globalPosY * outputL + globalPosX] + globalPosZ] = AddSample; - } - } + float l_s = (float)globalPosX; // the currently calculated line + + float delayMultiplicator = delayMultiplicatorRaw / s_i; + gDest[globalPosY * inputL + globalPosX] = delayMultiplicator * pow((l_s - l_i), 2) + s_i + (1-isPAImage)*s_i; + } +} - __kernel void ckDelayCalculationSphe( __global unsigned short *gDest, +__kernel void ckDelayCalculationSphe( __global unsigned short *gDest, __global unsigned short *usedLines, __global unsigned int *memoryLocations, unsigned int inputL, unsigned int inputS, unsigned int outputL, unsigned int outputS, char isPAImage, float delayMultiplicatorRaw // parameters ) - { - uint globalPosX = get_global_id(0); - uint globalPosY = get_global_id(1); - uint globalPosZ = get_global_id(2); - - if(globalPosZ < usedLines[globalPosY * 3 * outputL + 3 * globalPosX]) - { - float l_i = (float)globalPosX / outputL * inputL; - float s_i = (float)globalPosY / outputS * inputS / 2; - - float l_s = usedLines[globalPosY * 3 * outputL + 3 * globalPosX + 1] + globalPosZ; +{ + uint globalPosX = get_global_id(0); + uint globalPosY = get_global_id(1); + uint globalPosZ = get_global_id(2); + + if (globalPosX < inputL && globalPosY < outputS) + { + float l_i = 0; // we calculate the delays relative to line zero + float s_i = (float)globalPosY / (float)outputS * (float)inputS / 2; + + float l_s = (float)globalPosX; // the currently calculated line - gDest[memoryLocations[globalPosY * outputL + globalPosX] + globalPosZ] = - sqrt( - pow(s_i, 2) - + - pow((delayMultiplicatorRaw * ((l_s - l_i)) / inputL), 2) - ) + (1-isPAImage)*s_i; - } - } \ No newline at end of file + gDest[globalPosY * outputL + globalPosX] = + sqrt( + pow(s_i, 2) + + + pow((delayMultiplicatorRaw * ((l_s - l_i)) / inputL), 2) + ) + (1-isPAImage)*s_i; + } +} \ No newline at end of file diff --git a/Modules/PhotoacousticsAlgorithms/include/OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.h b/Modules/PhotoacousticsAlgorithms/include/OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.h index 38bee14c68..a947fe7583 100644 --- a/Modules/PhotoacousticsAlgorithms/include/OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.h +++ b/Modules/PhotoacousticsAlgorithms/include/OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.h @@ -1,105 +1,102 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #ifndef _MITKPHOTOACOUSTICSDELAYCALC_H_ #define _MITKPHOTOACOUSTICSDELAYCALC_H_ #ifdef PHOTOACOUSTICS_USE_GPU #include "mitkOclDataSetToDataSetFilter.h" #include #include "mitkPhotoacousticBeamformingSettings.h" namespace mitk { /** Documentation * * \brief The OclBinaryThresholdImageFilter computes a binary segmentation based on given threshold values. * * The filter requires two threshold values ( the upper and the lower threshold ) and two image values ( inside and outside ). The resulting voxel of the segmentation image is assigned the inside value 1 if the image value is between the given thresholds and the outside value otherwise. */ class OCLDelayCalculation : public OclDataSetToDataSetFilter, public itk::Object { public: mitkClassMacroItkParent(OCLDelayCalculation, itk::Object); itkNewMacro(Self); /** * @brief SetInput Set the input image. Only 3D images are supported for now. * @param image a 3D image. * @throw mitk::Exception if the dimesion is not 3. */ /** Update the filter */ void Update(); void SetConfig(BeamformingSettings conf) { m_Conf = conf; } - void SetInputs(cl_mem usedLines, cl_mem memoryLocations, unsigned int bufferSize) + void SetInputs(cl_mem usedLines) { - m_MemLoc = memoryLocations; m_UsedLines = usedLines; - m_BufferSize = bufferSize; } protected: /** Constructor */ OCLDelayCalculation(); /** Destructor */ virtual ~OCLDelayCalculation(); /** Initialize the filter */ bool Initialize(); void Execute(); mitk::PixelType GetOutputType() { return mitk::MakeScalarPixelType(); } int GetBytesPerElem() { return sizeof(unsigned short); } virtual us::Module* GetModule(); int m_sizeThis; private: /** The OpenCL kernel for the filter */ cl_kernel m_PixelCalculation; BeamformingSettings m_Conf; - cl_mem m_MemLoc; cl_mem m_UsedLines; unsigned int m_BufferSize; float m_DelayMultiplicatorRaw; char m_IsPAImage; }; } #endif #endif \ No newline at end of file diff --git a/Modules/PhotoacousticsAlgorithms/source/OpenCLFilter/mitkPhotoacousticOCLBeamformingFilter.cpp b/Modules/PhotoacousticsAlgorithms/source/OpenCLFilter/mitkPhotoacousticOCLBeamformingFilter.cpp index ef63c9e3e9..b529bea54a 100644 --- a/Modules/PhotoacousticsAlgorithms/source/OpenCLFilter/mitkPhotoacousticOCLBeamformingFilter.cpp +++ b/Modules/PhotoacousticsAlgorithms/source/OpenCLFilter/mitkPhotoacousticOCLBeamformingFilter.cpp @@ -1,266 +1,250 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #ifdef PHOTOACOUSTICS_USE_GPU #include "./OpenCLFilter/mitkPhotoacousticOCLBeamformingFilter.h" #include "usServiceReference.h" mitk::PhotoacousticOCLBeamformingFilter::PhotoacousticOCLBeamformingFilter() : m_PixelCalculation( NULL ), m_InputImage(mitk::Image::New()), m_ApodizationBuffer(nullptr), m_MemoryLocationsBuffer(nullptr), m_DelaysBuffer(nullptr), m_UsedLinesBuffer(nullptr) { this->AddSourceFile("DASQuadratic.cl"); this->AddSourceFile("DAS.cl"); this->AddSourceFile("DMAS.cl"); this->m_FilterID = "OpenCLBeamformingFilter"; this->Initialize(); unsigned int dim[] = { 128, 2048, 2 }; mitk::Vector3D spacing; spacing[0] = 1; spacing[1] = 1; spacing[2] = 1; m_InputImage->Initialize(mitk::MakeScalarPixelType(), 3, dim); m_InputImage->SetSpacing(spacing); m_ChunkSize[0] = 128; m_ChunkSize[1] = 128; m_ChunkSize[2] = 8; - m_SumFilter = mitk::OCLMemoryLocSum::New(); m_UsedLinesCalculation = mitk::OCLUsedLinesCalculation::New(); m_DelayCalculation = mitk::OCLDelayCalculation::New(); } mitk::PhotoacousticOCLBeamformingFilter::~PhotoacousticOCLBeamformingFilter() { if ( this->m_PixelCalculation ) { clReleaseKernel( m_PixelCalculation ); } if (m_ApodizationBuffer) clReleaseMemObject(m_ApodizationBuffer); } void mitk::PhotoacousticOCLBeamformingFilter::Update() { //Check if context & program available if (!this->Initialize()) { us::ServiceReference ref = GetModuleContext()->GetServiceReference(); OclResourceService* resources = GetModuleContext()->GetService(ref); // clean-up also the resources resources->InvalidateStorage(); mitkThrow() <<"Filter is not initialized. Cannot update."; } else{ // Execute this->Execute(); } } void mitk::PhotoacousticOCLBeamformingFilter::UpdateDataBuffers() { - if (BeamformingSettings::OutputDimensionsChanged(m_Conf,m_ConfOld)) + //Initialize the Output + try { - //Initialize the Output - try - { - MITK_INFO << "Updating Workgroup size for new dimensions"; - size_t outputSize = m_Conf.ReconstructionLines * m_Conf.SamplesPerLine * m_Conf.inputDim[2]; - m_OutputDim[0] = m_Conf.ReconstructionLines; - m_OutputDim[1] = m_Conf.SamplesPerLine; - m_OutputDim[2] = m_Conf.inputDim[2]; - this->InitExec(this->m_PixelCalculation, m_OutputDim, outputSize, sizeof(float)); - } - catch (const mitk::Exception& e) - { - MITK_ERROR << "Caught exception while initializing filter: " << e.what(); - return; - } + MITK_INFO << "Updating Workgroup size for new dimensions"; + size_t outputSize = m_Conf.ReconstructionLines * m_Conf.SamplesPerLine * m_Conf.inputDim[2]; + m_OutputDim[0] = m_Conf.ReconstructionLines; + m_OutputDim[1] = m_Conf.SamplesPerLine; + m_OutputDim[2] = m_Conf.inputDim[2]; + this->InitExec(this->m_PixelCalculation, m_OutputDim, outputSize, sizeof(float)); + } + catch (const mitk::Exception& e) + { + MITK_ERROR << "Caught exception while initializing filter: " << e.what(); + return; } if (BeamformingSettings::SettingsChangedOpenCL(m_Conf, m_ConfOld)) { cl_int clErr = 0; MITK_INFO << "Updating Buffers for new configuration"; // create the apodisation buffer if (m_Apodisation == nullptr) { MITK_INFO << "No apodisation function set; Beamforming will be done without any apodisation."; m_Apodisation = new float[1]; m_Apodisation[0] = 1; m_ApodArraySize = 1; } us::ServiceReference ref = GetModuleContext()->GetServiceReference(); OclResourceService* resources = GetModuleContext()->GetService(ref); cl_context gpuContext = resources->GetContext(); if (m_ApodizationBuffer) clReleaseMemObject(m_ApodizationBuffer); this->m_ApodizationBuffer = clCreateBuffer(gpuContext, CL_MEM_READ_ONLY | CL_MEM_USE_HOST_PTR, sizeof(float) * m_ApodArraySize, m_Apodisation, &clErr); CHECK_OCL_ERR(clErr); // calculate used lines m_UsedLinesCalculation->SetConfig(m_Conf); m_UsedLinesCalculation->Update(); - - // calculate memory locations of the Delays - - unsigned int sumDimensions[2] = { m_Conf.ReconstructionLines, m_Conf.SamplesPerLine }; - - m_SumFilter->SetInput(m_UsedLinesCalculation->GetGPUOutput()); - m_SumFilter->SetInputDimensions(sumDimensions); - m_SumFilter->Update(); - - unsigned int size = m_SumFilter->GetSum(); - m_UsedLinesBuffer = m_UsedLinesCalculation->GetGPUOutput()->GetGPUBuffer(); - m_MemoryLocationsBuffer = m_SumFilter->GetGPUOutput()->GetGPUBuffer(); // calculate the Delays m_DelayCalculation->SetConfig(m_Conf); - m_DelayCalculation->SetInputs(m_UsedLinesBuffer, m_MemoryLocationsBuffer, size); + m_DelayCalculation->SetInputs(m_UsedLinesBuffer); m_DelayCalculation->Update(); m_DelaysBuffer = m_DelayCalculation->GetGPUOutput()->GetGPUBuffer(); m_ConfOld = m_Conf; } } void mitk::PhotoacousticOCLBeamformingFilter::Execute() { cl_int clErr = 0; UpdateDataBuffers(); clErr = clSetKernelArg(this->m_PixelCalculation, 2, sizeof(cl_mem), &(this->m_UsedLinesBuffer)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 3, sizeof(cl_mem), &(this->m_MemoryLocationsBuffer)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 4, sizeof(cl_mem), &(this->m_DelaysBuffer)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 5, sizeof(cl_mem), &(this->m_ApodizationBuffer)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 6, sizeof(cl_ushort), &(this->m_ApodArraySize)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 7, sizeof(cl_uint), &(this->m_Conf.inputDim[0])); - clErr |= clSetKernelArg(this->m_PixelCalculation, 8, sizeof(cl_uint), &(this->m_Conf.inputDim[1])); - clErr |= clSetKernelArg(this->m_PixelCalculation, 9, sizeof(cl_uint), &(this->m_Conf.inputDim[2])); - clErr |= clSetKernelArg(this->m_PixelCalculation, 10, sizeof(cl_uint), &(this->m_Conf.ReconstructionLines)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 11, sizeof(cl_uint), &(this->m_Conf.SamplesPerLine)); - - CHECK_OCL_ERR(clErr); + clErr |= clSetKernelArg(this->m_PixelCalculation, 3, sizeof(cl_mem), &(this->m_DelaysBuffer)); + clErr |= clSetKernelArg(this->m_PixelCalculation, 4, sizeof(cl_mem), &(this->m_ApodizationBuffer)); + clErr |= clSetKernelArg(this->m_PixelCalculation, 5, sizeof(cl_ushort), &(this->m_ApodArraySize)); + clErr |= clSetKernelArg(this->m_PixelCalculation, 6, sizeof(cl_uint), &(this->m_Conf.inputDim[0])); + clErr |= clSetKernelArg(this->m_PixelCalculation, 7, sizeof(cl_uint), &(this->m_Conf.inputDim[1])); + clErr |= clSetKernelArg(this->m_PixelCalculation, 8, sizeof(cl_uint), &(this->m_Conf.inputDim[2])); + clErr |= clSetKernelArg(this->m_PixelCalculation, 9, sizeof(cl_uint), &(this->m_Conf.ReconstructionLines)); + clErr |= clSetKernelArg(this->m_PixelCalculation, 10, sizeof(cl_uint), &(this->m_Conf.SamplesPerLine)); // execute the filter on a 3D NDRange if (m_OutputDim[2] == 1 || m_ChunkSize[2] == 1) this->ExecuteKernelChunks(m_PixelCalculation, 2, m_ChunkSize); else this->ExecuteKernelChunks(m_PixelCalculation, 3, m_ChunkSize); // signalize the GPU-side data changed m_Output->Modified( GPU_DATA ); } us::Module *mitk::PhotoacousticOCLBeamformingFilter::GetModule() { return us::GetModuleContext()->GetModule(); } bool mitk::PhotoacousticOCLBeamformingFilter::Initialize() { bool buildErr = true; cl_int clErr = 0; if ( OclFilter::Initialize() ) { switch (m_Conf.Algorithm) { case BeamformingSettings::BeamformingAlgorithm::DAS: { this->m_PixelCalculation = clCreateKernel(this->m_ClProgram, "ckDAS", &clErr); break; } case BeamformingSettings::BeamformingAlgorithm::DMAS: { this->m_PixelCalculation = clCreateKernel(this->m_ClProgram, "ckDMAS", &clErr); break; } default: { MITK_INFO << "No beamforming algorithm specified, setting to DAS"; this->m_PixelCalculation = clCreateKernel(this->m_ClProgram, "ckDAS", &clErr); break; } } buildErr |= CHECK_OCL_ERR( clErr ); } + + CHECK_OCL_ERR(clErr); + return (OclFilter::IsInitialized() && buildErr ); } void mitk::PhotoacousticOCLBeamformingFilter::SetInput(mitk::Image::Pointer image) { OclDataSetToDataSetFilter::SetInput(image); m_InputImage = image; m_Conf.inputDim[0] = m_InputImage->GetDimension(0); m_Conf.inputDim[1] = m_InputImage->GetDimension(1); m_Conf.inputDim[2] = m_InputImage->GetDimension(2); } void mitk::PhotoacousticOCLBeamformingFilter::SetInput(void* data, unsigned int* dimensions, unsigned int BpE) { OclDataSetToDataSetFilter::SetInput(data, dimensions[0] * dimensions[1] * dimensions[2], BpE); m_Conf.inputDim[0] = dimensions[0]; m_Conf.inputDim[1] = dimensions[1]; m_Conf.inputDim[2] = dimensions[2]; } mitk::Image::Pointer mitk::PhotoacousticOCLBeamformingFilter::GetOutputAsImage() { mitk::Image::Pointer outputImage = mitk::Image::New(); if (m_Output->IsModified(GPU_DATA)) { void* pData = m_Output->TransferDataToCPU(m_CommandQue); const unsigned int dimension = 3; unsigned int dimensions[3] = { (unsigned int)m_OutputDim[0], (unsigned int)m_OutputDim[1], (unsigned int)m_OutputDim[2] }; const mitk::SlicedGeometry3D::Pointer p_slg = m_InputImage->GetSlicedGeometry(); MITK_DEBUG << "Creating new MITK Image."; outputImage->Initialize(this->GetOutputType(), dimension, dimensions); outputImage->SetSpacing(p_slg->GetSpacing()); outputImage->SetImportVolume(pData, 0, 0, mitk::Image::ReferenceMemory); } MITK_DEBUG << "Image Initialized."; return outputImage; } void* mitk::PhotoacousticOCLBeamformingFilter::GetOutput() { return OclDataSetToDataSetFilter::GetOutput(); } #endif diff --git a/Modules/PhotoacousticsAlgorithms/source/OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.cpp b/Modules/PhotoacousticsAlgorithms/source/OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.cpp index 13b216d8d8..32be8f2b86 100644 --- a/Modules/PhotoacousticsAlgorithms/source/OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.cpp +++ b/Modules/PhotoacousticsAlgorithms/source/OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.cpp @@ -1,119 +1,119 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ #define _USE_MATH_DEFINES #include #include "./OpenCLFilter/mitkPhotoacousticOCLDelayCalculation.h" #include "usServiceReference.h" #include "mitkImageReadAccessor.h" mitk::OCLDelayCalculation::OCLDelayCalculation() : m_PixelCalculation(NULL) { this->AddSourceFile("DelayCalculation.cl"); this->m_FilterID = "DelayCalculation"; this->Initialize(); } mitk::OCLDelayCalculation::~OCLDelayCalculation() { if (this->m_PixelCalculation) { clReleaseKernel(m_PixelCalculation); } } void mitk::OCLDelayCalculation::Update() { //Check if context & program available if (!this->Initialize()) { us::ServiceReference ref = GetModuleContext()->GetServiceReference(); OclResourceService* resources = GetModuleContext()->GetService(ref); // clean-up also the resources resources->InvalidateStorage(); mitkThrow() << "Filter is not initialized. Cannot update."; } else { // Execute this->Execute(); } } void mitk::OCLDelayCalculation::Execute() { cl_int clErr = 0; - unsigned int gridDim[3] = { m_Conf.ReconstructionLines, m_Conf.SamplesPerLine, m_Conf.ReconstructionLines }; + unsigned int gridDim[3] = { m_Conf.inputDim[0], m_Conf.SamplesPerLine, 1 }; + m_BufferSize = gridDim[0] * gridDim[1] * 1; try { this->InitExecNoInput(this->m_PixelCalculation, gridDim, m_BufferSize, sizeof(unsigned short)); } catch (const mitk::Exception& e) { MITK_ERROR << "Caught exception while initializing filter: " << e.what(); return; } if (m_Conf.DelayCalculationMethod == BeamformingSettings::DelayCalc::QuadApprox) m_DelayMultiplicatorRaw = pow(1 / (m_Conf.TimeSpacing*m_Conf.SpeedOfSound) * m_Conf.Pitch * m_Conf.TransducerElements / m_Conf.inputDim[0], 2) / 2; else if (m_Conf.DelayCalculationMethod == BeamformingSettings::DelayCalc::Spherical) m_DelayMultiplicatorRaw = 1 / (m_Conf.TimeSpacing*m_Conf.SpeedOfSound) * (m_Conf.Pitch*m_Conf.TransducerElements); m_IsPAImage = m_Conf.isPhotoacousticImage; - clErr |= clSetKernelArg(this->m_PixelCalculation, 1, sizeof(cl_mem), &m_UsedLines); - clErr |= clSetKernelArg(this->m_PixelCalculation, 2, sizeof(cl_mem), &m_MemLoc); - clErr |= clSetKernelArg(this->m_PixelCalculation, 3, sizeof(cl_uint), &(this->m_Conf.inputDim[0])); - clErr |= clSetKernelArg(this->m_PixelCalculation, 4, sizeof(cl_uint), &(this->m_Conf.inputDim[1])); - clErr |= clSetKernelArg(this->m_PixelCalculation, 5, sizeof(cl_uint), &(this->m_Conf.ReconstructionLines)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 6, sizeof(cl_uint), &(this->m_Conf.SamplesPerLine)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 7, sizeof(cl_char), &(this->m_IsPAImage)); - clErr |= clSetKernelArg(this->m_PixelCalculation, 8, sizeof(cl_float), &(this->m_DelayMultiplicatorRaw)); + clErr = clSetKernelArg(this->m_PixelCalculation, 1, sizeof(cl_mem), &(this->m_UsedLines)); + clErr |= clSetKernelArg(this->m_PixelCalculation, 2, sizeof(cl_uint), &(this->m_Conf.inputDim[0])); + clErr |= clSetKernelArg(this->m_PixelCalculation, 3, sizeof(cl_uint), &(this->m_Conf.inputDim[1])); + clErr |= clSetKernelArg(this->m_PixelCalculation, 4, sizeof(cl_uint), &(this->m_Conf.ReconstructionLines)); + clErr |= clSetKernelArg(this->m_PixelCalculation, 5, sizeof(cl_uint), &(this->m_Conf.SamplesPerLine)); + clErr |= clSetKernelArg(this->m_PixelCalculation, 6, sizeof(cl_char), &(this->m_IsPAImage)); + clErr |= clSetKernelArg(this->m_PixelCalculation, 7, sizeof(cl_float), &(this->m_DelayMultiplicatorRaw)); CHECK_OCL_ERR(clErr); // execute the filter on a 3D NDRange - this->ExecuteKernel(m_PixelCalculation, 3); + this->ExecuteKernel(m_PixelCalculation, 2); // signalize the GPU-side data changed m_Output->Modified(GPU_DATA); } us::Module *mitk::OCLDelayCalculation::GetModule() { return us::GetModuleContext()->GetModule(); } bool mitk::OCLDelayCalculation::Initialize() { bool buildErr = true; cl_int clErr = 0; if (OclFilter::Initialize()) { if(m_Conf.DelayCalculationMethod == BeamformingSettings::DelayCalc::QuadApprox) this->m_PixelCalculation = clCreateKernel(this->m_ClProgram, "ckDelayCalculationQuad", &clErr); if (m_Conf.DelayCalculationMethod == BeamformingSettings::DelayCalc::Spherical) this->m_PixelCalculation = clCreateKernel(this->m_ClProgram, "ckDelayCalculationSphe", &clErr); buildErr |= CHECK_OCL_ERR(clErr); } return (OclFilter::IsInitialized() && buildErr); } \ No newline at end of file