function  mask_ind = nic_icatb_createMask(files, fMRIHInfo, dataType, maskType, complexInfo)
%% Creates a mask based on the files given and mask type. If default mask is
% specified then first file of different subjects is loaded and indices greater than
% mean of all the voxels are taken into account. If mask is specified then
% only the indices in the mask are used.
%
% Input:
% 1. files - data structure with the field name
% 2. fMRIDIM - data dimensions
% 3. dataType - 'real' or 'complex'
% 4. complexInfo - structure containing file naming and complexType
%
% Output:
% 1. mask_ind - mask Indices



if ~exist('dataType', 'var')
    dataType = 'real';
end

if ~exist('maskType', 'var')
    maskType = 'default';
else
    if isempty(maskType)
        maskType = 'default';
    end
end

defaultMaskOption='first_file';

if isempty(defaultMaskOption)
    defaultMaskOption = 'first_file';
end


%% Get the extension (checks only the first data set)
% [pathstr, fName, extn] = fileparts(icatb_parseExtn(files(1).name(1, :)));
% 
% if (~strcmpi(extn, '.nii')  && ~strcmpi(extn, '.img'))
%     error('Mask should be specified in analyze or nifti format');
% end

disp('Creating Mask');

% Initialise data count
countN = 0;

%% Mask type
if strcmpi(maskType, 'default')
    if strcmpi(defaultMaskOption, 'all_files')
         disp('Using all files of the subjects to create default mask ...');
    else
         disp('Using first file of each subject to create default mask ...');
    end
else
    disp(['Using file ', deblank(files(1).name(1, :)), ' to create mask ...']);
end
complexInfo = struct('complex_file_naming',{'R_'  'I_'},'complexType','real&imaginary');

%% Default mask
if strcmpi(maskType, 'default')
        %% Loop over number of files
        for i = 1:length(files)
            
            %% Rename 4D nifti files
            tempF = nic_icatb_rename_4d_file(files{i}.name);

            %% Use first file or all files
            if ~strcmpi(defaultMaskOption, 'all_files')
                % use only the first time point
                tempF = deblank(tempF(1, :));
            end

            %% Include in-brain voxels
            for nn = 1:size(tempF, 1)
                countN = countN + 1;
                %% Load data
                data  = icatb_loadData(tempF, dataType, complexInfo, 'read', nn);
                if ~isreal(data)
                    data = abs(data);
                end

                %% Use voxels that surpass or equal mean
                data = data(:);
                tempInd = (data >= mean(data));

                if (countN == 1)
                    nonZeroInd = tempInd;
                else
                    nonZeroInd = (nonZeroInd) & tempInd;
                end

                clear tempInd data;
            end
            %% End for including in-brain voxels

            clear tempF;
            nic_waitbar(i/length(files),sprintf('Creating mask...'), 'Multiple-runs ICA (Preprocessing)','Child');
        end
        %% End for loop over number of files



else

    %% Calculate non-zero indices for selected mask
    dims = fMRIHInfo.V(1).dim(1:3); % Data dimensions
    data = nic_loadData(deblank(files{1}.name(1, :))); % Mask data
    size_data = size(data);
    if length(size_data) == 2
        size_data(3) = 1;
    end

    %% Check the dimensions of the mask w.r.t data
    if length(find(size_data == dims)) ~= length(dims)
        error('Error:MaskDim', 'Mask dimensions ([%s]) doesn''t match that of data dimensions ([%s])', ...
            num2str(size_data), num2str(dims));
    end

    %% Convert data to 2D matrix
    data = data(:);

    %% Replace Nan with zeros
    data(isnan(data)) = 0;
    nonZeroInd = (data ~= 0);
    clear data;

end
%% End for checking default or selected mask

%% Mask out the zeros
mask_ind = find(nonZeroInd ~= 0);

dims = fMRIHInfo.V(1).dim(1:3);
a=zeros(dims);
a(mask_ind)=1;

bwx=bwlabeln(a,26);
d=unique(bwx);
for i=1:length(d)-1
    ff=find(bwx==i);
    if(length(ff)<500)
        a(ff)=0;
    end
end
mask_ind = find(a ~= 0);

if isempty(mask_ind)
    error('Error in creating Mask');
end

disp('Done Creating Mask');

return;


function files = nic_icatb_rename_4d_file(files, fileNumber)
%For nifti and analyze files rename the files by adding a number at the
%end

if ~exist('fileNumber', 'var')
    fileNumber = [];
end

if ~isempty(files)
    
    % Do pattern match to check IMG or NII files
    files = cellstr(files);
    
    % MATCH IMG OR NII
    checkNII = regexpi(files, '\.nii$');
    checkIMG = regexpi(files, '\.img$');
    
    good_inds1 = nic_icatb_good_cells(checkNII);
    good_inds2 = nic_icatb_good_cells(checkIMG);
    
    % Get good cells
    good_inds = (good_inds1 | good_inds2);
    good_inds = find(good_inds ~= 0);
    
    clear good_inds1 good_inds2 checkNII checkIMG;
    
    % If IMG or NII files exist check the headers of the image files
    if ~isempty(good_inds)
        imFiles = files(good_inds);
        if (~isempty(fileNumber)) && (length(fileNumber) == 1) && (fileNumber == 1)
            imFiles = strcat(imFiles, ',1');
            files(good_inds) = imFiles;
        else
            % Loop over img or nii files
            for nIm = 1:length(imFiles)
                currentFile = imFiles{nIm};
                
                try
                    ni = read_nifti(currentFile);
                    numFiles = ni.dat.dim(4);
                catch
                    numFiles = 1;
                end
                
                if (isempty(fileNumber))
                    tempFileNum = (1:numFiles);
                else
                    tempFileNum = fileNumber;
                end
                
                tempFileNum(tempFileNum > numFiles) = [];
                
                if ~isempty(tempFileNum)
                    tempFiles = [repmat(currentFile, length(tempFileNum), 1), repmat(',', length(tempFileNum), 1), num2str(tempFileNum(:))];
                else
                    tempFiles = '';
                end
                files{good_inds(nIm)} = tempFiles;
            end
            % End loop over img or nii files
            files = cellstr(str2mat(files));
            ind = nic_icatb_good_cells(files);
            files = files(ind);
        end
    end
    % End for checking IMG or NII files
    
    if ~isempty(files)
        files = str2mat(files);
    end
    
end

function ind = nic_icatb_good_cells(mycell)
%% Find good cells

if ~iscell(mycell)
    mycell = {mycell};
end

ind = cellfun('isempty', mycell);

% Good cells
ind = (ind == 0);

% Convert ind to row vector
ind = ind(:)';
