nonhuman primate take 1

src/MATLAB/+RSF/darkTubes.m

@@ -8,6 +8,7 @@ vx = max(voxelSize)./voxelSize;
 sigma = 1/sqrt(3) * (thickness ./ 2) * vx;
 
 cudaTarget = getCudaTarget;
+im = mat2gray(im);
 imf = HIP.LoG(im,sigma,cudaTarget);
 
 if thicknessIn < 0
@@ -15,11 +16,13 @@ if thicknessIn < 0
     imn = imf;
     imn(imn>0)=0;
     imn=abs(imn);
-    sigmaPlate = [5,1,1].*sigma;
+    sigmaPlate = [1,10,1].*sigma;
     imfPlate = HIP.LoG(imn,sigmaPlate,cudaTarget);
     imnPlate = imfPlate;
     imnPlate(imnPlate>0) = 0;
     imnPlate = abs(imnPlate);
+    imnPlate = imnPlate .* imbinarize(imnPlate);
+%     imnPlate = mat2gray(imnPlate);
     imnPlate = imnPlate .* imregionalmax(imnPlate);
     imOut = imnPlate;
 else        

src/MATLAB/nonHumanPrimate/goHuman.m

+tic
+path(path,fullfile('../clusterExamples/'))
+% ROOT = '/g/leverjs/Schuman_OCT/OCT/nonHumanPrimate';
+% fname = fullfile(ROOT,'31344_OS_V_5x5_0_0004257.LEVER');
+% ROOT = '/g/leverjs/Schuman_OCT/OCT/nonHumanPrimate';
+% fname = fullfile(ROOT,'31344_OS_V_5x5_0_0004257.LEVER');
+
+
+
+[imraw,CONSTANTS] = leversc.loadImage(fname,1,2);
+% imraw = flip(permute(imraw,[2,1]),2); % ACK verify direction!
+rx = 5;
+im_medfilt = medfilt3(imraw,[5,5,5]);
+% note NLM takes minutes (vs. e.g. seconds for LoG etc)...use batch!
+cudaTarget = getCudaTarget;
+im = HIP.NLMeans(im_medfilt,0.1,12,1,cudaTarget);
+
+thickness = 5;
+
+voxelSize = CONSTANTS.imageData.PixelPhysicalSize;
+voxelSize([1,2]) = voxelSize([2,1]); % swap (x,y,z) to matlab column major (y,x,z)
+vx = max(voxelSize)./voxelSize;
+sigma = 1/sqrt(3) * thickness./2 * vx;
+
+cudaTarget = getCudaTarget;
+
+imf = HIP.LoG(im,sigma,cudaTarget);
+imp=imf;imp(imp<0)=0;
+
+imn = imf;
+imn(imn>0)=0;imn=abs(imn);
+
+sigmaPlate = [1,5,1].*sigma;
+imfPlate = HIP.LoG(imn,sigmaPlate,cudaTarget);
+imnPlate = imfPlate;
+imnPlate(imnPlate>0) = 0;
+imnPlate = abs(imnPlate);
+imnPlate = mat2gray(imnPlate);
+
+voxelSize = CONSTANTS.imageData.PixelPhysicalSize;
+voxelSize([1,2]) = voxelSize([2,1]); % swap (x,y,z) to matlab column major (y,x,z)
+vx = max(voxelSize)./voxelSize;
+sigma = 1/sqrt(3) * thickness./2 * vx;
+
+cudaTarget = getCudaTarget;
+imf = HIP.LoG(im,sigma,cudaTarget);
+imn = imf;
+imn(imn>0) = 0;
+imn = abs(imn);
+imp=imf;imp(imp<0)=0;
+imfm = fibermetric((imp));
+
+outfile = 'nhp_metricEmbeddings.tif';
+
+%
+imfmX = imfm.*imregionalmax(imfm);
+
+
+
+% bright plate
+sigmaPlate = [1,5,1].*sigma;
+imfPlate = HIP.LoG(im,sigmaPlate,cudaTarget);
+
+imnPlate = imfPlate;
+imnPlate(imnPlate>0) = 0;
+imnPlate = abs(imnPlate);
+imnPlateX = imnPlate.*imregionalmax(imnPlate,6);
+imnPlateX = mat2gray(imnPlateX);
+% impPlate = imfPlate;
+% impPlate(imnPlate<0) = 0;
+% impPlateX = impPlate.*imregionalmax(imnPlate,6);
+
+sliceThickness = 10;
+bClear = true;
+for Z = 100: 1   :size(im,1)-100
+    
+    clipLimits=[];
+    clipLimits.dim = 1;
+    clipLimits.range = [Z-sliceThickness:Z+sliceThickness];
+    clipLimits.range(clipLimits.range>size(im,1))=[];
+    clipLimits.range(clipLimits.range<1)=[];
+    
+    im1 = getClipFrame(im,clipLimits,gray());
+    im2 = getClipFrame(imn,clipLimits);
+    im3 = getClipFrame(imp,clipLimits);
+    im4 = getClipFrame(imfm,clipLimits);    
+    im6 = getClipFrame(imraw,clipLimits,gray());
+%     im7 = getClipFrame(imfmX,clipLimits);
+    
+    im7 = getClipFrame(imnPlateX,clipLimits);
+    im7 = imdilate(im7,strel('disk',3));
+    im8 = getClipFrame(imnPlate,clipLimits);
+
+    imOut = [ im6,im2, im8 ; im3, im4, im7   ];
+    clf;imagesc(imOut);
+    title(num2str(Z))
+    drawnow
+    if bClear
+        imwrite(imOut,outfile);
+        bClear = false;
+    else
+        imwrite(imOut,outfile,'writemode','append');
+    end
+end
+toc
\ No newline at end of file

src/MATLAB/nonHumanPrimate/goNHP.m

+tic
+path(path,fullfile('../clusterExamples/'))
+ROOT = '/g/leverjs/Schuman_OCT/OCT/nonHumanPrimate';
+fname = fullfile(ROOT,'31344_OS_V_5x5_0_0004257.LEVER');
+
+[imraw,CONSTANTS] = leversc.loadImage(fname,1,2);
+% imraw = flip(permute(imraw,[2,1]),2); % ACK verify direction!
+rx = 5;
+im_medfilt = medfilt3(imraw,[3,5,5]); % experimental denoise 
+% note NLM takes minutes (vs. e.g. seconds for LoG etc)...use batch!
+cudaTarget = getCudaTarget;
+im = HIP.NLMeans(im_medfilt,0.1,12,1,cudaTarget);
+
+thickness = 5;
+
+voxelSize = CONSTANTS.imageData.PixelPhysicalSize;
+voxelSize([1,2]) = voxelSize([2,1]); % swap (x,y,z) to matlab column major (y,x,z)
+vx = max(voxelSize)./voxelSize;
+sigma = 1/sqrt(3) * thickness./2 * vx;
+
+cudaTarget = getCudaTarget;
+
+imf = HIP.LoG(im,sigma,cudaTarget);
+imp=imf;imp(imp<0)=0;
+
+imn = imf;
+imn(imn>0)=0;imn=abs(imn);
+
+sigmaPlate = [1,5,1].*sigma;
+imfPlate = HIP.LoG(imn,sigmaPlate,cudaTarget);
+imnPlate = imfPlate;
+imnPlate(imnPlate>0) = 0;
+imnPlate = abs(imnPlate);
+imnPlate = mat2gray(imnPlate);
+
+voxelSize = CONSTANTS.imageData.PixelPhysicalSize;
+voxelSize([1,2]) = voxelSize([2,1]); % swap (x,y,z) to matlab column major (y,x,z)
+vx = max(voxelSize)./voxelSize;
+sigma = 1/sqrt(3) * thickness./2 * vx;
+
+cudaTarget = getCudaTarget;
+imf = HIP.LoG(im,sigma,cudaTarget);
+imn = imf;
+imn(imn>0) = 0;
+imn = abs(imn);
+imp=imf;imp(imp<0)=0;
+imfm = fibermetric((imp));
+
+outfile = 'nhp_metricEmbeddings.tif';
+
+%
+imfmX = imfm.*imregionalmax(imfm);
+
+
+
+% bright plate
+sigmaPlate = [1,5,1].*sigma;
+imfPlate = HIP.LoG(im,sigmaPlate,cudaTarget);
+
+imnPlate = imfPlate;
+imnPlate(imnPlate>0) = 0;
+imnPlate = abs(imnPlate);
+imnPlateX = imnPlate.*imregionalmax(imnPlate,6);
+imnPlateX = mat2gray(imnPlateX);
+% impPlate = imfPlate;
+% impPlate(imnPlate<0) = 0;
+% impPlateX = impPlate.*imregionalmax(imnPlate,6);
+
+sliceThickness = 10;
+bClear = true;
+for Z = 100: 1   :size(im,1)-100
+    
+    clipLimits=[];
+    clipLimits.dim = 1;
+    clipLimits.range = [Z-sliceThickness:Z+sliceThickness];
+    clipLimits.range(clipLimits.range>size(im,1))=[];
+    clipLimits.range(clipLimits.range<1)=[];
+    
+    im1 = getClipFrame(im,clipLimits,gray());
+    im2 = getClipFrame(imn,clipLimits);
+    im3 = getClipFrame(imp,clipLimits);
+    im4 = getClipFrame(imfm,clipLimits);    
+    im6 = getClipFrame(imraw,clipLimits,gray());
+%     im7 = getClipFrame(imfmX,clipLimits);
+    
+    im7 = getClipFrame(imnPlateX,clipLimits);
+    im7 = imdilate(im7,strel('disk',3));
+    im8 = getClipFrame(imnPlate,clipLimits);
+
+    imOut = [ im6,im2, im8 ; im3, im4, im7   ];
+    clf;imagesc(imOut);
+    title(num2str(Z))
+    drawnow
+    if bClear
+        imwrite(imOut,outfile);
+        bClear = false;
+    else
+        imwrite(imOut,outfile,'writemode','append');
+    end
+end
+toc
\ No newline at end of file