Commit 30c085bc authored by Eric Wait's avatar Eric Wait

Added a segmentor

parent d13b79ee
......@@ -294,6 +294,15 @@ for i=1:7
fprintf('%s took %f sec and returned a dif of %f\n',kernelName,sumTime,dif);
clear imageOut;
tic
kernelName = 'Segment';
imageOut = CudaMex(sprintf('%s',kernalName),image1,1.0,[5,5,5],device);
fprintf('%s took %f sec\n',kernelName,toc);
if (showOut)
showIm(imageOut,[kernelName ' Max']);
end
clear imageOut;
threshold = additive;
tic
kernelName = 'ThresholdFilter';
......
......@@ -15,6 +15,7 @@
#include "CudaPow.cuh"
#include "CudaImageReduction.cuh"
#include "CudaSum.cuh"
#include "CudaSegment.cuh"
#include "CudaThreshold.cuh"
void clearDevice()
......@@ -814,6 +815,42 @@ void clearDevice()
}
unsigned char* segment(const unsigned char* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel/*=NULL*/, unsigned char** imageOut/*=NULL*/, int device/*=0*/)
{
return cSegment(imageIn,dims,alpha,kernelDims,kernel,imageOut,device);
}
unsigned short* segment(const unsigned short* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel/*=NULL*/, unsigned short** imageOut/*=NULL*/, int device/*=0*/)
{
return cSegment(imageIn,dims,alpha,kernelDims,kernel,imageOut,device);
}
short* segment(const short* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel/*=NULL*/, short** imageOut/*=NULL*/, int device/*=0*/)
{
return cSegment(imageIn,dims,alpha,kernelDims,kernel,imageOut,device);
}
unsigned int* segment(const unsigned int* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel/*=NULL*/, unsigned int** imageOut/*=NULL*/, int device/*=0*/)
{
return cSegment(imageIn,dims,alpha,kernelDims,kernel,imageOut,device);
}
int* segment(const int* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel/*=NULL*/, int** imageOut/*=NULL*/, int device/*=0*/)
{
return cSegment(imageIn,dims,alpha,kernelDims,kernel,imageOut,device);
}
float* segment(const float* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel/*=NULL*/, float** imageOut/*=NULL*/, int device/*=0*/)
{
return cSegment(imageIn,dims,alpha,kernelDims,kernel,imageOut,device);
}
double* segment(const double* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel/*=NULL*/, double** imageOut/*=NULL*/, int device/*=0*/)
{
return cSegment(imageIn,dims,alpha,kernelDims,kernel,imageOut,device);
}
unsigned char* thresholdFilter(const unsigned char* imageIn, Vec<size_t> dims, unsigned char thresh, unsigned char** imageOut/*=NULL*/, int device/*=0*/)
{
return cThresholdFilter(imageIn,dims,thresh,imageOut,device);
......@@ -848,3 +885,4 @@ void clearDevice()
{
return cThresholdFilter(imageIn,dims,thresh,imageOut,device);
}
......@@ -190,9 +190,17 @@ IMAGE_PROCESSOR_API size_t sumArray(const int* imageIn, size_t n, int device=0);
IMAGE_PROCESSOR_API double sumArray(const float* imageIn, size_t n, int device=0);
IMAGE_PROCESSOR_API double sumArray(const double* imageIn, size_t n, int device=0);
IMAGE_PROCESSOR_API unsigned char* segment(const unsigned char* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel=NULL, unsigned char** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API unsigned short* segment(const unsigned short* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel=NULL, unsigned short** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API short* segment(const short* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel=NULL, short** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API unsigned int* segment(const unsigned int* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel=NULL, unsigned int** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API int* segment(const int* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel=NULL, int** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API float* segment(const float* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel=NULL, float** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API double* segment(const double* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel=NULL, double** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API unsigned char* thresholdFilter(const unsigned char* imageIn, Vec<size_t> dims, unsigned char thresh, unsigned char** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API unsigned short* thresholdFilter(const unsigned short* imageIn, Vec<size_t> dims, unsigned short thresh, unsigned short** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API short* thresholdFilter(const short* imageIn, Vec<size_t> dims, int thresh, short** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API short* thresholdFilter(const short* imageIn, Vec<size_t> dims, short thresh, short** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API unsigned int* thresholdFilter(const unsigned int* imageIn, Vec<size_t> dims, unsigned int thresh, unsigned int** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API int* thresholdFilter(const int* imageIn, Vec<size_t> dims, int thresh, int** imageOut=NULL, int device=0);
IMAGE_PROCESSOR_API float* thresholdFilter(const float* imageIn, Vec<size_t> dims, float thresh, float** imageOut=NULL, int device=0);
......
#pragma once
#include "CudaThreshold.cuh"
#include "CudaMaxFilter.cuh"
#include "CudaMinFilter.cuh"
template <class PixelType>
PixelType* cSegment(const PixelType* imageIn, Vec<size_t> dims, double alpha, Vec<size_t> kernelDims, float* kernel=NULL,
PixelType** imageOut=NULL, int device=0)
{
PixelType thresh = cOtsuThresholdValue(imageIn,dims,device);
thresh = (PixelType)(thresh*alpha);
PixelType* imOut = setUpOutIm(dims, imageOut);
PixelType minVal = std::numeric_limits<PixelType>::min();
PixelType maxVal = std::numeric_limits<PixelType>::max();
cudaDeviceProp props;
cudaGetDeviceProperties(&props,device);
if (kernel==NULL)
{
kernelDims = kernelDims.clamp(Vec<size_t>(1,1,1),dims);
float* ones = new float[kernelDims.product()];
for (int i=0; i<kernelDims.product(); ++i)
ones[i] = 1.0f;
HANDLE_ERROR(cudaMemcpyToSymbol(cudaConstKernel, ones, sizeof(float)*kernelDims.product()));
delete[] ones;
}
else
{
HANDLE_ERROR(cudaMemcpyToSymbol(cudaConstKernel, kernel, sizeof(float)*kernelDims.product()));
}
size_t availMem, total;
cudaMemGetInfo(&availMem,&total);
std::vector<ImageChunk> chunks = calculateBuffers<PixelType>(dims,2,(size_t)(availMem*MAX_MEM_AVAIL),props,kernelDims);
Vec<size_t> maxDeviceDims;
setMaxDeviceDims(chunks, maxDeviceDims);
CudaDeviceImages<PixelType> deviceImages(2,maxDeviceDims,device);
DEBUG_KERNEL_CHECK();
for (std::vector<ImageChunk>::iterator curChunk=chunks.begin(); curChunk!=chunks.end(); ++curChunk)
{
curChunk->sendROI(imageIn,dims,deviceImages.getCurBuffer());
deviceImages.setNextDims(curChunk->getFullChunkSize());
cudaThreshold<<<curChunk->blocks,curChunk->threads>>>(*(deviceImages.getCurBuffer()),*(deviceImages.getNextBuffer()),thresh,minVal,maxVal);
DEBUG_KERNEL_CHECK();
deviceImages.incrementBuffer();
cudaMinFilter<<<curChunk->blocks,curChunk->threads>>>(*(deviceImages.getCurBuffer()),*(deviceImages.getNextBuffer()),kernelDims,
minVal,maxVal);
DEBUG_KERNEL_CHECK();
deviceImages.incrementBuffer();
cudaMaxFilter<<<curChunk->blocks,curChunk->threads>>>(*(deviceImages.getCurBuffer()),*(deviceImages.getNextBuffer()),kernelDims,
minVal,maxVal);
DEBUG_KERNEL_CHECK();
deviceImages.incrementBuffer();
curChunk->retriveROI(imOut,dims,deviceImages.getCurBuffer());
}
return imOut;
}
\ No newline at end of file
......@@ -42,6 +42,7 @@
<ClInclude Include="Cuda\CudaNormalizedCovariance.cuh" />
<ClInclude Include="Cuda\CudaPolyTransferFunc.cuh" />
<ClInclude Include="Cuda\CudaPow.cuh" />
<ClInclude Include="Cuda\CudaSegment.cuh" />
<ClInclude Include="Cuda\CudaSum.cuh" />
<ClInclude Include="Cuda\CudaThreshold.cuh" />
<ClInclude Include="Cuda\CudaUtilities.cuh" />
......
......@@ -104,6 +104,9 @@
<ClInclude Include="Cuda\ImageChunk.cuh">
<Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="Cuda\CudaSegment.cuh">
<Filter>Header Files</Filter>
</ClInclude>
</ItemGroup>
<ItemGroup>
<ClCompile Include="Cuda\CHelpers.cpp">
......
......@@ -42,6 +42,7 @@
<ClInclude Include="Cuda\CudaNormalizedCovariance.cuh" />
<ClInclude Include="Cuda\CudaPolyTransferFunc.cuh" />
<ClInclude Include="Cuda\CudaPow.cuh" />
<ClInclude Include="Cuda\CudaSegment.cuh" />
<ClInclude Include="Cuda\CudaSum.cuh" />
<ClInclude Include="Cuda\CudaThreshold.cuh" />
<ClInclude Include="Cuda\CudaUtilities.cuh" />
......
......@@ -104,6 +104,9 @@
<ClInclude Include="Cuda\ImageChunk.cuh">
<Filter>Header Files</Filter>
</ClInclude>
<ClInclude Include="Cuda\CudaSegment.cuh">
<Filter>Header Files</Filter>
</ClInclude>
</ItemGroup>
<ItemGroup>
<ClCompile Include="Cuda\CHelpers.cpp">
......
......@@ -127,6 +127,7 @@ copy $(OutDir)\CudaMex.dll ..\..\src\MATLAB\CudaMex.mexw64</Command>
<ClCompile Include="Mex\MexOtsuThresholdFilter.cpp" />
<ClCompile Include="Mex\MexOtsuThresholdValue.cpp" />
<ClCompile Include="Mex\MexReduceImage.cpp" />
<ClCompile Include="Mex\MexSegment.cpp" />
<ClCompile Include="Mex\MexSumArray.cpp" />
<ClCompile Include="Mex\MexThresholdFilter.cpp" />
</ItemGroup>
......
......@@ -113,6 +113,9 @@
<ClCompile Include="Mex\MexThresholdFilter.cpp">
<Filter>Source Files</Filter>
</ClCompile>
<ClCompile Include="Mex\MexSegment.cpp">
<Filter>Source Files</Filter>
</ClCompile>
</ItemGroup>
<ItemGroup>
<None Include="Mex\CudaMex.def">
......
......@@ -35,6 +35,7 @@ void MexCommand::init()
REGISTER_COMMAND(MexOtsuThresholdValue);
REGISTER_COMMAND(MexReduceImage);
REGISTER_COMMAND(MexSumArray);
REGISTER_COMMAND(MexSegment);
REGISTER_COMMAND(MexThresholdFilter);
}
......
......@@ -335,9 +335,22 @@ public:
virtual void execute(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]);
virtual std::string check(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]);
virtual std::string printUsage();
virtual std::string printHelp();
virtual std::string printHelp();
};
class MexSegment : MexCommand
{
public:
MexSegment(){}
virtual ~MexSegment(){}
virtual void execute(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]);
virtual std::string check(int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[]);
virtual std::string printUsage();
virtual std::string printHelp();
};
class MexReduceImage : MexCommand
{
public:
......
#include "MexCommand.h"
#include "Vec.h"
#include "CWrappers.cuh"
#include "CHelpers.h"
void MexSegment::execute( int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[] )
{
int device = 0;
if (nrhs>3)
device = mat_to_c((int)mxGetScalar(prhs[3]));
double alpha = mxGetScalar(prhs[1]);
double* radiusMex = (double*)mxGetData(prhs[2]);
Vec<size_t> radius((int)radiusMex[0],(int)radiusMex[1],(int)radiusMex[2]);
Vec<size_t> kernDims;
float* kern = createEllipsoidKernel(radius,kernDims);
Vec<size_t> imageDims;
if (mxIsUint8(prhs[0]))
{
unsigned char* imageIn,* imageOut;
setupImagePointers(prhs[0],&imageIn,&imageDims,&plhs[0],&imageOut);
segment(imageIn,imageDims,alpha,kernDims,kern,&imageOut,device);
}
else if (mxIsUint16(prhs[0]))
{
unsigned short* imageIn,* imageOut;
setupImagePointers(prhs[0],&imageIn,&imageDims,&plhs[0],&imageOut);
segment(imageIn,imageDims,alpha,kernDims,kern,&imageOut,device);
}
else if (mxIsInt16(prhs[0]))
{
short* imageIn,* imageOut;
setupImagePointers(prhs[0],&imageIn,&imageDims,&plhs[0],&imageOut);
segment(imageIn,imageDims,alpha,kernDims,kern,&imageOut,device);
}
else if (mxIsUint32(prhs[0]))
{
unsigned int* imageIn,* imageOut;
setupImagePointers(prhs[0],&imageIn,&imageDims,&plhs[0],&imageOut);
segment(imageIn,imageDims,alpha,kernDims,kern,&imageOut,device);
}
else if (mxIsInt32(prhs[0]))
{
int* imageIn,* imageOut;
setupImagePointers(prhs[0],&imageIn,&imageDims,&plhs[0],&imageOut);
segment(imageIn,imageDims,alpha,kernDims,kern,&imageOut,device);
}
else if (mxIsSingle(prhs[0]))
{
float* imageIn,* imageOut;
setupImagePointers(prhs[0],&imageIn,&imageDims,&plhs[0],&imageOut);
segment(imageIn,imageDims,alpha,kernDims,kern,&imageOut,device);
}
else if (mxIsDouble(prhs[0]))
{
double* imageIn,* imageOut;
setupImagePointers(prhs[0],&imageIn,&imageDims,&plhs[0],&imageOut);
segment(imageIn,imageDims,alpha,kernDims,kern,&imageOut,device);
}
else
{
mexErrMsgTxt("Image type not supported!");
}
}
std::string MexSegment::check( int nlhs, mxArray* plhs[], int nrhs, const mxArray* prhs[] )
{
if (nrhs<3 || nrhs>4)
return "Incorrect number of inputs!";
if (nlhs!=1)
return "Requires one output!";
size_t numDims = mxGetNumberOfDimensions(prhs[0]);
if (numDims>3)
return "Image can have a maximum of three dimensions!";
if (!mxIsDouble(prhs[1]))
return "Alpha has to be a single double!";
size_t numEl = mxGetNumberOfElements(prhs[2]);
if (numEl!=3 || !mxIsDouble(prhs[2]))
return "Median neighborhood has to be an array of three doubles!";
return "";
}
std::string MexSegment::printUsage()
{
return "imageOut = CudaMex('Segment',imageIn,alpha,[MorphClosureX,MorphClosureY,MorphClosureZ],[device]);";
}
std::string MexSegment::printHelp()
{
std::string msg = "\tSegmentaion is done by applying an Otsu adaptive threshold (which can be modified by the alpha multiplier).\n";
msg += "\tA morphological closing is then applied using a ellipsoid neighborhood with the MorphClosure dimensions.\n";
msg += "\n";
return msg;
}
\ No newline at end of file
Markdown is supported
0%
or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment