diff --git a/fslmaths.cc b/fslmaths.cc
index 0dc7c3f8023bf068c52cda27da681293f76e7d8d..d5b6a163ebb4aea0644d9625e076b42702f22277 100755
--- a/fslmaths.cc
+++ b/fslmaths.cc
@@ -11,7 +11,7 @@
#include "newimage/newimageall.h"
#include "miscmaths/miscmaths.h"
#include "utils/fsl_isfinite.h"
-#include "libprob/libprob.h"
+#include "libprob.h"
using namespace MISCMATHS;
using namespace NEWIMAGE;
@@ -65,6 +65,7 @@ int printUsage(const string& programName)
cout << " -bin : use (current image>0) to binarise" << endl;
cout << " -binv : binarise and invert (binarisation and logical inversion)" << endl;
cout << " -fillh : fill holes in a binary mask (holes are internal - i.e. do not touch the edge of the FOV)" << endl;
+ cout << " -fillh26 : fill holes using 26 connectivity" << endl;
cout << " -index : replace each nonzero voxel with a unique (subject to wrapping) index number" << endl;
cout << " -grid <value> <spacing> : add a 3D grid of intensity <value> with grid spacing <spacing>" << endl;
cout << " -edge : edge strength" << endl;
@@ -84,8 +85,9 @@ int printUsage(const string& programName)
cout << "\nKernel operations (set BEFORE filtering operation if desired):" << endl;
cout << " -kernel 3D : 3x3x3 box centered on target voxel (set as default kernel)" << endl;
cout << " -kernel 2D : 3x3x1 box centered on target voxel" << endl;
- cout << " -kernel box <size> : all voxels in a box of width <size> centered on target voxel" << endl;
- cout << " -kernel boxv <size> : <size>x<size>x<size> box centered on target voxel, CAUTION: size should be an odd number" << endl;
+ cout << " -kernel box <size> : all voxels in a cube of width <size> mm centered on target voxel" << endl;
+ cout << " -kernel boxv <size> : all voxels in a cube of width <size> voxels centered on target voxel, CAUTION: size should be an odd number" << endl;
+ cout << " -kernel boxv3 <X> <Y> <Z>: all voxels in a cuboid of dimensions X x Y x Z centered on target voxel, CAUTION: size should be an odd number" << endl;
cout << " -kernel gauss <sigma> : gaussian kernel (sigma in mm, not voxels)" << endl;
cout << " -kernel sphere <size> : all voxels in a sphere of radius <size> mm centered on target voxel" << endl;
cout << " -kernel file <filename> : use external file as kernel" << endl;
@@ -149,7 +151,7 @@ void loadNewImage(volume4D<T> &oldI, volume4D<T> &newI, string filename)
volume4D<T> tmpvol=oldI;
oldI.reinitialize(oldI.xsize(),oldI.ysize(),oldI.zsize(),newI.tsize());
oldI.copyproperties(tmpvol);
- for (int t=0;t<newI.tsize();t++) oldI[t]=tmpvol.constantSubVolume(0);
+ for (int t=0;t<newI.tsize();t++) oldI[t]=tmpvol[0];
}
// sanity check on valid orientation info (it should be consistent)
@@ -184,6 +186,8 @@ int inputParser(int argc, char *argv[], short output_dt)
read_volume4D(inputVolume,string(argv[1]));
bool modifiedInput(false);
bool setDisplayRange(false);
+ float tfce_delta(0);
+ float tfce_minT(0);
int i=2;
for (i = 2; i < argc-1; i++) //main loop
@@ -216,7 +220,12 @@ int inputParser(int argc, char *argv[], short output_dt)
else if (string(argv[i]+2) == "min") inputVolume.value(x,y,z,t)=column_volume.min();
else if (string(argv[i]+2) == "mean") inputVolume.value(x,y,z,t)=(T)column_volume.mean();
else if (string(argv[i]+2) == "std") inputVolume.value(x,y,z,t)=(T)column_volume.stddev();
- else if (string(argv[i]+2) == "maxn") inputVolume.value(x,y,z,t)=column_volume.maxcoordx();
+ else if (string(argv[i]+2) == "maxn") {
+ vector<long long> iCoords;
+ volume<T> mask;
+ column_volume.max(mask,iCoords);
+ inputVolume.value(x,y,z,t)=iCoords[0];
+ }
else if (string(argv[i]+2) == "median") inputVolume.value(x,y,z,t)=column_volume.percentile(0.5);
else if (string(argv[i]+2) == "perc") inputVolume.value(x,y,z,t)=column_volume.percentile(atof(argv[i+1])/100.0);
else if (string(argv[i]+2) == "ar1") {
@@ -320,8 +329,7 @@ int truecount=(int)truth.sum();
// TPim
volume4D<float> TPim;
copyconvert(loginput,TPim);
-for(int t=0;t<loginput.tsize();t++)
- TPim[t]*=truth;
+TPim*=truth;
// noise
if (!separatenoise)
@@ -332,8 +340,7 @@ if (!separatenoise)
invtruth.value(x,y,z)=0;
invtruecount=(int)invtruth.sum();
copyconvert(loginput,noise);
- for(int t=0;t<loginput.tsize();t++)
- noise[t]*=invtruth;
+ noise*=invtruth;
}
// }}}
@@ -404,9 +411,8 @@ if (!separatenoise)
else if (string(argv[i])=="-mas")
{
loadNewImage(inputVolume, temp_volume, string(argv[++i]));
- temp_volume.binarise(0,temp_volume.max()+1,exclusive); // needed to binarise max() value + 1 due to
- for (int t=0;t<inputVolume.tsize();t++) //potential issue with exclusive binarise
- inputVolume[t]*=temp_volume[t%temp_volume.tsize()]; //this gives compatibility with 3 and 4D masks
+ temp_volume.binarise(0,temp_volume.max()+1,exclusive); // needed to binarise max() value + 1 due to potential issue with exclusive binarise
+ inputVolume*=temp_volume; //this gives compatibility with 3 and 4D masks
} //without needing to have a specific volume3D variable
/***************************************************************/
else if (string(argv[i])=="-add"){
@@ -415,7 +421,7 @@ if (!separatenoise)
else
{
loadNewImage(inputVolume, temp_volume, string(argv[i]));
- for (int t=0;t<inputVolume.tsize();t++) inputVolume[t]+=temp_volume[t%temp_volume.tsize()];
+ inputVolume+=temp_volume;
}}
/***************************************************************/
else if (string(argv[i])=="-sub"){
@@ -424,7 +430,7 @@ if (!separatenoise)
else
{
loadNewImage(inputVolume, temp_volume, string(argv[i]));
- for (int t=0;t<inputVolume.tsize();t++) inputVolume[t]-=temp_volume[t%temp_volume.tsize()];
+ inputVolume-=temp_volume;
}}
/***************************************************************/
else if (string(argv[i])=="-mul"){
@@ -433,7 +439,7 @@ if (!separatenoise)
else
{
loadNewImage(inputVolume, temp_volume, string(argv[i]));
- for (int t=0;t<inputVolume.tsize();t++) inputVolume[t]*=temp_volume[t%temp_volume.tsize()];
+ inputVolume*=temp_volume;
}}
/***************************************************************/
else if (string(argv[i])=="-rem"){
@@ -565,14 +571,16 @@ if (!separatenoise)
else
{
float size=atof(argv[i+2]);
- if(string(argv[i+1])=="box") kernel=box_kernel(size,xdim,ydim,zdim);
- if(string(argv[i+1])=="boxv") kernel=box_kernel((int)size,(int)size,(int)size);
+ separable=false;
+ if(string(argv[i+1])=="box" || string(argv[i+1])=="boxv" || string(argv[i+1])=="boxv3" || string(argv[i+1])=="gauss") separable=true;
+
+ if(string(argv[i+1])=="box") kernel=box_kernel(size,xdim,ydim,zdim);
+ else if(string(argv[i+1])=="boxv") kernel=box_kernel((int)size,(int)size,(int)size);
+ else if(string(argv[i+1])=="boxv3") { kernel=box_kernel((int)atof(argv[i+2]),(int)atof(argv[i+3]),(int)atof(argv[i+4])); i+=2; }
else if(string(argv[i+1])=="gauss") kernel=gaussian_kernel3D(size,xdim,ydim,zdim);
else if(string(argv[i+1])=="sphere") kernel=spherical_kernel(size,xdim,ydim,zdim);
else if(string(argv[i+1])=="file") read_volume(kernel,string(argv[i+2]));
- if(string(argv[i+1])=="box" || string(argv[i+1])=="gauss") separable=true;
- else separable=false;
- i++;
+ i++;
}
i++;
//save_volume(kernel,"kernel");
@@ -894,7 +902,16 @@ if (!separatenoise)
{
inputVolume.binarise(0,inputVolume.max()+1,exclusive);
for (int t=0; t<=inputVolume.maxt(); t++) {
- inputVolume[t] = fill_holes(inputVolume[t]);
+ inputVolume[t] = fill_holes(inputVolume[t],6);
+ }
+ inputVolume.binarise(0,inputVolume.max()+1,exclusive);
+ }
+ /***************************************************************/
+ else if (string(argv[i])=="-fillh26")
+ {
+ inputVolume.binarise(0,inputVolume.max()+1,exclusive);
+ for (int t=0; t<=inputVolume.maxt(); t++) {
+ inputVolume[t] = fill_holes(inputVolume[t],26);
}
inputVolume.binarise(0,inputVolume.max()+1,exclusive);
}
@@ -927,10 +944,11 @@ if (!separatenoise)
/*****************SPATIAL FILTERING OPTIONS*********************/
/***********************All Dilation***************************/
else if (string(argv[i])=="-dilall") {
- volume4D<T> mask(inputVolume);
- mask.binarise(0,0,inclusive); // get a zero mask
- mask=(T)1-mask; // convert to a non-zero mask
- for(int t=0;t<inputVolume.tsize();t++) dilall(inputVolume[t],mask[t]);
+ for(int t=0;t<inputVolume.tsize();t++) {
+ volume<T> mask(inputVolume[t]);
+ mask.binarise(0,0,inclusive,true); // get a zero mask, convert to non-zero
+ inputVolume[t]=dilall(inputVolume[t],mask);
+ }
}
/***********************Mean Dilation***************************/
else if (string(argv[i])=="-dilM")
@@ -959,6 +977,10 @@ if (!separatenoise)
/*****************END OF FILTERING OPTIONS***************/
else if (string(argv[i])=="-edge")
inputVolume=edge_strengthen(inputVolume);
+ else if (string(argv[i])=="-tfce_minT")
+ tfce_minT = atof(argv[++i]);
+ else if (string(argv[i])=="-tfce_delta")
+ tfce_delta = atof(argv[++i]);
else if (string(argv[i])=="-tfce") {
float height_power = atof(argv[++i]);
float size_power = atof(argv[++i]);
@@ -966,7 +988,9 @@ if (!separatenoise)
for(int t=0;t<inputVolume.tsize();t++) {
try {
- tfce(inputVolume[t], height_power, size_power, connectivity, 0, 0);
+ volume<T> im(inputVolume[t]);
+ tfce(im, height_power, size_power, connectivity, tfce_minT, tfce_delta);
+ inputVolume[t]=im;
}
catch(Exception& e) {
cerr << "ERROR: TFCE failed, please check your file for valid sizes and voxel values." << e.what() << endl << endl << "Exiting" << endl;
@@ -987,8 +1011,11 @@ if (!separatenoise)
int Z = atoi(argv[++i]);
float tfce_thresh = atof(argv[++i]);
- for(int t=0;t<inputVolume.tsize();t++)
- tfce_support(inputVolume[t], height_power, size_power, connectivity, 0, 0, X, Y, Z, tfce_thresh);
+ for(int t=0;t<inputVolume.tsize();t++) {
+ volume<T> im(inputVolume[t]);
+ tfce_support(im, height_power, size_power, connectivity, tfce_minT, tfce_delta, X, Y, Z, tfce_thresh);
+ inputVolume[t]=im;
+ }
}
// }}}
@@ -1127,7 +1154,7 @@ if (!separatenoise)
dti_V3(inputVolume.xsize(),inputVolume.ysize(),inputVolume.zsize(),3);
dti_L1=0; dti_L2=0; dti_L3=0; dti_FA=0; dti_MD=0; dti_MO=0; dti_V1=0; dti_V2=0; dti_V3=0;
-
+
copybasicproperties(inputVolume,dti_L1);
copybasicproperties(inputVolume,dti_L2);
copybasicproperties(inputVolume,dti_L3);
@@ -1197,14 +1224,23 @@ if (!separatenoise)
// if i+1>argc-1 then can save (otherwise it is a syntax error with no specific output specified)
check_for_output_name(i+1,argc-1);
+ dti_L1.setDisplayMaximumMinimum(dti_L1.max(),dti_L1.min());
save_volume(dti_L1,string(argv[argc-1])+"_L1");
+ dti_L2.setDisplayMaximumMinimum(dti_L2.max(),dti_L2.min());
save_volume(dti_L2,string(argv[argc-1])+"_L2");
+ dti_L3.setDisplayMaximumMinimum(dti_L3.max(),dti_L3.min());
save_volume(dti_L3,string(argv[argc-1])+"_L3");
+ dti_FA.setDisplayMaximumMinimum(1,0);
save_volume(dti_FA,string(argv[argc-1])+"_FA");
+ dti_MD.setDisplayMaximumMinimum(dti_MD.max(),dti_MD.min());
save_volume(dti_MD,string(argv[argc-1])+"_MD");
+ dti_MO.setDisplayMaximumMinimum(1,-1);
save_volume(dti_MO,string(argv[argc-1])+"_MO");
+ dti_V1.setDisplayMaximumMinimum(1,-1);
save_volume4D(dti_V1,string(argv[argc-1])+"_V1");
+ dti_V2.setDisplayMaximumMinimum(1,-1);
save_volume4D(dti_V2,string(argv[argc-1])+"_V2");
+ dti_V3.setDisplayMaximumMinimum(1,-1);
save_volume4D(dti_V3,string(argv[argc-1])+"_V3");
// }}}
@@ -1236,7 +1272,7 @@ if (!separatenoise)
inputVolume.setDisplayMaximumMinimum((float)inputVolume.max(),(float)inputVolume.min());
- save_volume4D_datatype(inputVolume,string(argv[argc-1]),output_dt);
+ save_volume_datatype(inputVolume,string(argv[argc-1]),output_dt);
return 0;
}