Merge pull request #23 from adhusch/develop

Regular Merge

lib/ANTSTransforms/applyANTSTransformToPoints.m

@@ -22,7 +22,15 @@
 % (c) 2017
 % [email protected], [email protected]
 function transPoints = applyANTSTransformToPoints(points, antsTransformFileStrings)
-applyTransformsToPointsCmd = '/Applications/ANTSGit/bin/antsApplyTransformsToPoints'; % FIXME: Make multi-os compatabile and find bins automatically
+try
+[~, t] = system('echo $ANTSPATH');
+antspath = t(1:end-1); % remove line ending
+end
+if(~isempty(t))
+    applyTransformsToPointsCmd = [antspath 'antsApplyTransformsToPoints']; 
+else
+    applyTransformsToPointsCmd = 'antsApplyTransformsToPoints'; % FIXME: Make multi-os compatabile and find bins automatically
+end
 
 tempFile = [tempname() '.csv'];
 tempFileOut = [tempname() '.csv'];

lib/NiftiOOP/NiftiMod.m

@@ -85,8 +85,10 @@
         function [image, header] = load(this)
             if( ~this.isLoaded && ~strcmp(this.filepath,'')) % it is not already isLoaded and we were not called with empty string (by Matlab callingn our get functions without command after initiating!)
                 disp(['Loading ' this.filepath ' from disk...']);
-                % nifti = load_nii(this.filepath, [], [], [], [], [], 0.1, 'q'); % load with 0.1 tolerance, prefer the q-form over the sform!
-                nifti = load_untouch_nii(this.filepath); % load with 0.1 tolerance, prefer the q-form over the sform!
+                %nifti = load_nii(this.filepath, [], [], [], [], [], 0.1, 'q'); % load with 0.1 tolerance, prefer the q-form over the sform!
+                nifti = load_untouch_nii(this.filepath); % load untouched, i.e. keep orginial data orientation (dont apply any transform)
+                nifti = applyNiiIntensityScaling(nifti); %however, the voxel intesities should be scaled correctly now here (as they
+%would by load_nii otherewise)
 
                 this.header  = nifti.hdr; % header should maybe kept all the time ..
 

lib/NiftiOOP/applyNiiIntensityScaling.m

@@ -0,0 +1,61 @@
+function nii = applyNiiIntensityScaling(nii)
+% Code extracte from Jimmy Shen's ([email protected]) xform_nii.m
+% method, which is part of the nifti toolbox
+%
+% Andreas Husch, 2019
+
+    %  if scl_slope field is nonzero, then each voxel value in the
+       %  dataset should be scaled as: y = scl_slope * x + scl_inter
+       %  I bring it here because hdr will be modified by change_hdr.
+       %
+       if nii.hdr.dime.scl_slope ~= 0 & ...
+        ismember(nii.hdr.dime.datatype, [2,4,8,16,64,256,512,768]) & ...
+        (nii.hdr.dime.scl_slope ~= 1 | nii.hdr.dime.scl_inter ~= 0)
+
+          nii.img = ...
+        nii.hdr.dime.scl_slope * double(nii.img) + nii.hdr.dime.scl_inter;
+
+          if nii.hdr.dime.datatype == 64
+
+             nii.hdr.dime.datatype = 64;
+             nii.hdr.dime.bitpix = 64;
+          else
+             nii.img = single(nii.img);
+
+             nii.hdr.dime.datatype = 16;
+             nii.hdr.dime.bitpix = 32;
+          end
+
+          nii.hdr.dime.glmax = max(double(nii.img(:)));
+          nii.hdr.dime.glmin = min(double(nii.img(:)));
+
+          %  set scale to non-use, because it is applied in xform_nii
+          %
+          nii.hdr.dime.scl_slope = 0;
+
+       end
+
+       %  However, the scaling is to be ignored if datatype is DT_RGB24.
+
+       %  If datatype is a complex type, then the scaling is to be applied
+       %  to both the real and imaginary parts.
+       %
+       if nii.hdr.dime.scl_slope ~= 0 & ...
+        ismember(nii.hdr.dime.datatype, [32,1792])
+
+          nii.img = ...
+        nii.hdr.dime.scl_slope * double(nii.img) + nii.hdr.dime.scl_inter;
+
+          if nii.hdr.dime.datatype == 32
+             nii.img = single(nii.img);
+          end
+
+          nii.hdr.dime.glmax = max(double(nii.img(:)));
+          nii.hdr.dime.glmin = min(double(nii.img(:)));
+
+          %  set scale to non-use, because it is applied in xform_nii
+          %
+          nii.hdr.dime.scl_slope = 0;
+
+       end
+end

src/Functions/extractElectrodePointclouds.m

@@ -50,9 +50,8 @@
     disp(['Thresholding ' niiCT.filepath  ' for metal with METAL_THRESHOLD = ' num2str(METAL_THRESHOLD) '...']);
     maskedImg = niiCT.img;
 
-   %  [xx,yy,zz] = ndgrid(-10:10);
-    %structEle = sqrt(xx.^2 + yy.^2 + zz.^2) <= 2.5 / sqrt(max(niiCT.voxsize));
-    %brainMask = imerode(imerode(brainMask,structEle),structEle);
+    structEle = strel('sphere', ceil(3 / max(niiCT.voxsize)) ); % make sure brain mask contains no skull
+    brainMask = imerode(brainMask,structEle);
 
     maskedImg(~(brainMask)) = NaN;
     threImg = (maskedImg > METAL_THRESHOLD);