combineROI

 combineROI -- combine two or three ROIs by intersection, union, or overlap

 combineROI takes a list of two or three ROIs and combines them by one
 of three methods:

   intersection: set intersection; if a voxel is present in all of the ROIs,
   the voxel is in the output set

   union: set union; if a voxel is present in one of the ROIs, it is present
   in the output set.

   overlap: voxels are tagged with a color according to whether or not 
   they overlap with voxels in another ROI.

 The proper way to call combineROI is as follows:

 ROI = combineROI('-method', 'intersection|union|overlap', ...
 'path_to_your_roi', {cell_list_of_up_to_three_rois}, ...
 'stem_of_output_roi', [number_of_slices number_of_rows number_of_columns]);

 The output variable "ROI" is the result of the combination.  Note that this
 is of size [number_of_slices number_of_rows number_of_columns] and commonly
 should be compacted using compactROI.
 
 $Id: combineROI.m,v 1.8 2003/09/15 15:52:18 nknouf Exp $

0001 function [ROI] = combineROI(varargin)
0002 % combineROI -- combine two or three ROIs by intersection, union, or overlap
0003 %
0004 % combineROI takes a list of two or three ROIs and combines them by one
0005 % of three methods:
0006 %
0007 %   intersection: set intersection; if a voxel is present in all of the ROIs,
0008 %   the voxel is in the output set
0009 %
0010 %   union: set union; if a voxel is present in one of the ROIs, it is present
0011 %   in the output set.
0012 %
0013 %   overlap: voxels are tagged with a color according to whether or not
0014 %   they overlap with voxels in another ROI.
0015 %
0016 % The proper way to call combineROI is as follows:
0017 %
0018 % ROI = combineROI('-method', 'intersection|union|overlap', ...
0019 % 'path_to_your_roi', {cell_list_of_up_to_three_rois}, ...
0020 % 'stem_of_output_roi', [number_of_slices number_of_rows number_of_columns]);
0021 %
0022 % The output variable "ROI" is the result of the combination.  Note that this
0023 % is of size [number_of_slices number_of_rows number_of_columns] and commonly
0024 % should be compacted using compactROI.
0025 %
0026 % $Id: combineROI.m,v 1.8 2003/09/15 15:52:18 nknouf Exp $
0027 
0028 % License:    Perl Artistic License
0029 % History:    9/10/03: initial release
0030 
0031 %% if there are no arguments, print usage
0032 if (nargin == 0)
0033     printUsage;
0034     return;
0035 end
0036 
0037 %% if there is only one argument, it means we should initialize the structure
0038 %% or show the help information
0039 if (nargin == 1)
0040     flag = deblank(varargin{1});
0041     if (strncmp(flag, '-help', 5))
0042         showHelp;
0043         return;
0044     else
0045         printUsage;
0046         return;
0047     end
0048 %% if there are six arguments, it means that we should
0049 %% actually combine the ROIs
0050 elseif (nargin == 6)
0051     flag = deblank(varargin{1});
0052     if (strncmp(flag, '-method', 7))
0053         %% if the second argument is "intersection", we compute
0054         %% the intersection
0055         flag = deblank(varargin{2});
0056         if (strncmp(flag, 'intersection', 12))
0057             ROI = createCompositeROI(varargin{3}, varargin{4}, 1, varargin{6});
0058         %% else, we compute the union
0059         elseif (strncmp(flag, 'union', 5))
0060             ROI = createCompositeROI(varargin{3}, varargin{4}, 2, varargin{6});
0061         elseif (strncmp(flag, 'overlap', 7))
0062             outputStem = deblank(varargin{5});
0063             infoFile = sprintf('%s%s.info', varargin{3}, outputStem);
0064             ROI = createOverlapMap(varargin{3}, varargin{4}, infoFile, varargin{6});
0065         else
0066             printUsage;
0067             return;
0068         end
0069 
0070         %% now let's save the returned ROI
0071         outputStem = deblank(varargin{5});
0072         stem = sprintf('%s%s.roi', varargin{3}, outputStem);
0073         roimtx = compactROI(ROI);
0074         if (isempty(roimtx))
0075             roimtx = [1 1 1 0];
0076         end
0077         s = sprintf('save %s roimtx -ASCII',stem);
0078         eval(s);
0079         return;
0080     else
0081         printUsage;
0082         return;
0083     end
0084 else
0085     printUsage;
0086 end
0087 
0088 
0089 %% this function combines the ROIs together using intersections or unions
0090 function ROI = createCompositeROI(roiPath, roiList, method, roiSize)
0091     %% get the number of ROIs
0092     [junk, numROIs] = size(roiList);
0093     if (numROIs == 1)
0094         oops(252);
0095     end
0096 
0097     %% get the dimensions of the ROIs
0098     stem = sprintf('%s%s', roiPath, roiList{1,1});
0099     %%[nrows, ncols, ntp, fs, ns, endian, bext] = fmri_bfiledim(stem);
0100     ns = roiSize(1,1);
0101     nrows = roiSize(1,2);
0102     ncols = roiSize(1,3);
0103 
0104     stem = sprintf('%s%s.roi', roiPath, roiList{1,1});
0105 
0106     %% create empty output ROI
0107     ROI = zeros(ns, nrows, ncols);
0108 
0109 
0110     %% let's premptively load the first ROI in roiList
0111     %% this will simplify things later
0112     ROI = textread(stem);
0113     ROI = expandROI(ROI, [ns nrows ncols]);
0114     %% change to column major format
0115     ROI = permute(ROI, [1 3 2]);
0116     
0117     %% we need to do pairwise comparisons, so we loop from
0118     %% 1 to numROIs-1
0119     for i=1:(numROIs-1)
0120         %% now, we need to load the i+1th volume
0121         stem = sprintf('%s%s.roi', roiPath, roiList{1,i+1});
0122         vol = textread(stem);
0123         vol = expandROI(vol, [ns nrows ncols]);
0124         %% change to column major format
0125         vol = permute(vol, [1 3 2]);
0126         %%[vol] = read_bfloat(stem, ns);
0127 
0128         %% next, either do the intersection or union
0129         if (method == 1)
0130             ROI = ROI & vol;    
0131         elseif (method == 2)
0132             ROI = ROI | vol;
0133         end    
0134     end
0135 return;
0136 
0137 %% create a map of the overlap of ROIs
0138 function ROI = createOverlapMap(roiPath, roiList, infoFile, roiSize)
0139     %% get the number of ROIs
0140     [junk, numROIs] = size(roiList);
0141     if (numROIs == 1)
0142         oops(252);
0143     elseif (numROIs > 3)
0144         oops(253);
0145     end
0146 
0147     %% get the dimensions of the ROIs
0148     stem = sprintf('%s%s', roiPath, roiList{1,1});
0149     %%[nrows, ncols, ntp, fs, ns, endian, bext] = fmri_bfiledim(stem);
0150     ns = roiSize(1,1);
0151     nrows = roiSize(1,2);
0152     ncols = roiSize(1,3);
0153 
0154     %% create empty output ROI
0155     ROI = zeros(ns, nrows, ncols);
0156 
0157     %% load in all of the ROIs
0158     %% while we're at it, why don't we save the
0159     %% list of the ROIs used in this overlap.
0160     f = fopen(infoFile, 'w');
0161     if (~f)
0162         oops(158);
0163     end
0164     for i=1:numROIs
0165         stem = sprintf('%s%s.roi', roiPath, roiList{1,i});
0166         ROIs{i} = textread(stem);
0167         ROIs{i} = expandROI(ROIs{i}, [ns nrows ncols]); 
0168         %% change to column major format
0169         ROIs{i} = permute(ROIs{i}, [1 3 2]);
0170         %%ROIs{i} = read_bfloat(stem, ns);
0171         fprintf(f,'%s\n', roiList{1,i});
0172     end
0173     fclose(f);
0174 
0175     fprintf('\nComputing the overlap of the ROIs\n');
0176 
0177     %% big loop; should find more effecient way to write this later
0178     for i=1:ns
0179         fprintf('On slice %d\n', i);
0180     for j=1:nrows
0181     for k=1:ncols
0182         %% setup empty temp variable
0183         temp = [];
0184 
0185         %% create vector of the values in the ROIs
0186         for m=1:numROIs
0187             newtemp = ROIs{m};
0188             temp = [temp squeeze(newtemp(i,j,k))];
0189         end
0190     
0191         %% only do the following when we have non-zero values
0192         %% in at least one of the ROIs
0193         while (find(temp>0))
0194             switch numROIs
0195             %% need to check which ROIs were in the overlap
0196             %% and color-code the voxels appropriately.
0197             %% This is trivial in the case of 2 ROIs, more
0198             %% complicated in the case of 3, and 4?  Well,
0199             %% we won't even go there.  If anyone wants to
0200             %% implement more than 3 ROIs, be my guest.
0201             case 2,
0202                 if ((temp(1,1) == 1) & (temp(1,2) == 1))
0203                     ROI(i,j,k) = 4;
0204                 elseif (temp(1,1) == 1)
0205                     ROI(i,j,k) = 2;
0206                 elseif (temp(1,2) == 1)
0207                     ROI(i,j,k) = 3;
0208                 end
0209                 break;
0210             case 3,
0211                 if ((temp(1,1) == temp(1,2)) & (temp(1,1) == temp(1,3)))
0212                     ROI(i,j,k) = 8;
0213                 elseif ((temp(1,1) == 1) & (temp(1,2) == 1))
0214                     ROI(i,j,k) = 5;
0215                 elseif ((temp(1,1) == 1) & (temp(1,3) == 1))
0216                     ROI(i,j,k) = 6;
0217                 elseif ((temp(1,2) == 1) & (temp(1,3) == 1))
0218                     ROI(i,j,k) = 7;
0219                 elseif (temp(1,1) == 1)
0220                     ROI(i,j,k) = 2;
0221                 elseif (temp(1,2) == 1)
0222                     ROI(i,j,k) = 3;
0223                 elseif (temp(1,3) == 1)
0224                     ROI(i,j,k) = 4;
0225                 end
0226                 break;
0227             end
0228         end
0229     end
0230     end
0231     end
0232 return;
0233 
0234 %% print basic usage information when no parameters are given
0235 function printUsage()
0236     fprintf(1, '%%%%%%%%%%  combineROI.m  %%%%%%%%%%\n');
0237     fprintf(1, ['$Id: combineROI.m,v 1.8 2003/09/15 15:52:18 nknouf Exp $' '\n\n']);
0238     fprintf(1, 'computeROI(''-help'') -- gives detailed help information\n');
0239     fprintf(1, 'ROI = combineROI(''-method'', ''<intersection | union>'', roiPath, [roiList], outputStem) -- combines the ROIs in ROIlist by either intersection or union\n');
0240 
0241     return;
0242 
0243 function showHelp()
0244     fprintf(1, '$Id: combineROI.m,v 1.8 2003/09/15 15:52:18 nknouf Exp $\n\n');
0245     fprintf(1, 'this will show some help information someday\n');
0246 
0247     return;