%
%
%  Time Projection Function for Imaris
%
%  Copyright Bitplane AG 2006
%
%
%  Installation:
%
%  - Copy this file into the XTensions folder in the Imaris installation directory.
%  - You will find this function in the Image Processing menu
%
%    <CustomTools>
%      <Menu>
%        <Item name="Time Projection" icon="Matlab" tooltip="Project all channels of a time dataset to a single time point.">
%          <Command>Matlab::BPTimeProjection(%i)</Command>
%        </Item>
%      </Menu>
%    </CustomTools>
%  
%
%  Description:
%   
%   Project all channels of a time dataset to a single time point. 
%   For each voxel the intensity is computed using a MIP or mean projection
%       of all time points.
%   The projection mode is asked to the user. 
%
%

function BPTimeProjection(aImarisApplicationID, aProjectionMode)

% vImarisApplication = actxserver('Imaris.Application'); vImarisApplication.mVisible = true;

if isa(aImarisApplicationID, 'COM.Imaris_Application')
    vImarisApplication = aImarisApplicationID;
else
    % connect to Imaris Com interface
    vImarisServer = actxserver('ImarisServer.Server');
    vImarisApplication = vImarisServer.GetObject(aImarisApplicationID);
end
vImarisApplication.DataSetPushUndo('Time Projection');

% get the data set 
vImarisDataSet = vImarisApplication.mDataSet.Clone;

vNumberOfTimePoints = vImarisDataSet.mSizeT;
vNumberOfChannels = vImarisDataSet.mSizeC;
vNumberOfSlices = vImarisDataSet.mSizeZ;

% MIP or Mean projection
if nargin<2
    vAnswer = questdlg('Please select the projection mode:', 'Time Projection', ...
        'MIP', 'Mean','Cancel', 'MIP');
    if strcmp(vAnswer, 'Cancel')
        return;
    else
        vProjectionMode = vAnswer;
    end
else
    vProjectionMode = aProjectionMode;
end

vTotalCount = vNumberOfChannels*vNumberOfSlices;
vProgressDisplay = waitbar(0, 'Projecting...');
vProgressCount = 0;

vImage = zeros(vImarisDataSet.mSizeX, vImarisDataSet.mSizeY, vNumberOfTimePoints);

for vChannel = 1:vNumberOfChannels
    for vSlice = 1:vNumberOfSlices
        
        for vTime = 1:vNumberOfTimePoints
            vImage(:,:,vTime) = vImarisDataSet.GetDataSlice(vSlice-1, vChannel-1, vTime-1);
        end
        
        if strcmp(vProjectionMode, 'MIP')
            % MIP
            vImageSlice = max(vImage, [], 3); 
        else
            % Mean 
            vImageSlice = mean(vImage, 3);  
        end
        
        if strcmp(vImarisDataSet.mType,'eTypeUInt8')
            vImageSlice = uint8(vImageSlice);
        elseif strcmp(vImarisDataSet.mType,'eTypeUInt16')
            vImageSlice = uint16(vImageSlice);
        else
            vImageSlice = single(vImageSlice);
        end
        
        vImarisDataSet.SetDataSlice(vImageSlice, vSlice-1, vChannel-1, 0);
        vProgressCount = vProgressCount + 1;

        waitbar(vProgressCount/vTotalCount)
    end
end

vImarisApplication.mVisibleIndexT = 0;
vImarisDataSet.mSizeT = 1;

close(vProgressDisplay);

vImarisApplication.mDataSet = vImarisDataSet;