%
%
%  Track Plot Distance Between Tracks 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>
%       <Submenu name="Tracks Functions">
%        <Item name="Plot Distance Between Tracks" icon="Matlab" tooltip="Plot distance between tracks.">
%          <Command>Matlab::BPTrackPlotDistanceBetweenTracks(%i)</Command>
%        </Item>
%       </Submenu>
%      </Menu>
%      <SurpassTab>
%        <SurpassComponent name="bpTrack">
%          <Item name="Plot Distance Between Tracks" icon="Matlab" tooltip="Plot distance between tracks.">
%            <Command>Matlab::BPTrackPlotDistanceBetweenTracks(%i)</Command>
%          </Item>
%        </SurpassComponent>
%      </SurpassTab>
%    </CustomTools>
% 
%
%  Description:
%   
%   Plot the absolute value and the X,Y,Z components of the distance 
%       between the selected track and the other visible tracks that are 
%       in the same group.
%
%

function BPTrackPlotDistanceBetweenTracks(aImarisApplicationID)

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

% the user has to create a scene with some tracks
vSurpassScene = vImarisApplication.mSurpassScene;
if isequal(vSurpassScene, [])
    msgbox('Please create some tracks in the surpass scene!');
    return;
end

% get the track
vTrack1 = vImarisApplication.mFactory.ToTrack(vImarisApplication.mSurpassSelection);

if ~vImarisApplication.mFactory.IsTrack(vTrack1)  
    msgbox('Please select a track!');
    return;
end

% get the spots and edges
vSpots1 = vTrack1.GetSpots;
[vSpots1XYZ,vSpots1Time,vSpots1Radius] = vSpots1.Get;
vTimeMin1 = min(vSpots1Time);
vTimeMax1 = max(vSpots1Time);

vTimeRange1 = zeros(vTimeMax1-vTimeMin1+1,1);
vTimeRange1(vSpots1Time-vTimeMin1+1) = 1:length(vSpots1Time);

% get the second track
vParentGroup = vTrack1.GetParent;

vTrackNames = {};

% get the jet colormap
vColorMap = jet;
vColorMapSize = size(vColorMap,1);

for vPlot = 1:4
    subplot(2,2,vPlot)
    newplot
    hold on
end

vNumberOfTracks = 0;
vTracksList = [];

vNumberOfChildren = vParentGroup.GetNumberOfChildren;
for vChildIndex = 1:vNumberOfChildren
    vDataItem = vParentGroup.GetChild(vChildIndex - 1);

    % search the other visible tracks
    if vImarisApplication.mFactory.IsTrack(vDataItem) && ... 
            ~strcmp(vDataItem.mName,vTrack1.mName) && ... 
            vDataItem.mVisible == true
        
        vNumberOfTracks = vNumberOfTracks + 1;
        vTracksList(vNumberOfTracks) = vChildIndex;
        
    end
end

vNumberOfTracksDisplayed = 0;

for vTrackIndex = 1:vNumberOfTracks
    vTrack2 = vParentGroup.GetChild(vTracksList(vTrackIndex) - 1);

    vSpots2 = vTrack2.GetSpots;
    [vSpots2XYZ,vSpots2Time,vSpots2Radius] = vSpots2.Get;
    vTimeMin2 = min(vSpots2Time);
    vTimeMax2 = max(vSpots2Time);

    % both tracks are defined in [vTimeMin, vTimeMax]
    vTimeMin = max([vTimeMin1,vTimeMin2]);
    vTimeMax = min([vTimeMax1,vTimeMax2]);

    if vTimeMin>vTimeMax
        % the two tracks are not defined on a common time interval
        continue;
    end

    vSpotsTime = [];
    vSpots1Index = [];
    vSpots2Index = [];
    vNumberOfTimes = 0;
    
    % find the time points where both tracks are defined
    for vTime = vSpots1Time
        if vTime>=vTimeMin && vTime<=vTimeMax
            vTime2 = find(vSpots2Time==vTime);
            if length(vTime2)>0
                vNumberOfTimes = vNumberOfTimes + 1;
                vSpotsTime(vNumberOfTimes) = vTime + 1;
                vSpots1Index(vNumberOfTimes) = vTimeRange1(vTime-vTimeMin1+1);
                vSpots2Index(vNumberOfTimes) = vTime2(1);
            end
        end
    end
    
    if vNumberOfTimes<1
        % the two tracks are not commonly defined on any time point
        continue;
    end

    vDistanceX = vSpots1XYZ(vSpots1Index,1)-vSpots2XYZ(vSpots2Index,1);
    vDistanceY = vSpots1XYZ(vSpots1Index,2)-vSpots2XYZ(vSpots2Index,2);
    vDistanceZ = vSpots1XYZ(vSpots1Index,3)-vSpots2XYZ(vSpots2Index,3);
    vDistance = sqrt(vDistanceX.^2+vDistanceY.^2+vDistanceZ.^2);

    % computes the color
    vPercent = (vTrackIndex-1)/max(vNumberOfTracks-1,1)*(vColorMapSize-1)+1;
    vFloor = floor(vPercent);
    vColor = vColorMap(vFloor,:);
    if vFloor<vPercent
        vPercent = vPercent - vFloor;
        vColor = vPercent*vColor+(1-vPercent)*vColorMap(vFloor+1,:);
    end

    vStrings = {'','X','Y','Z'};
    for vPlot = 1:4
        subplot(2,2,vPlot)
        plot(vSpotsTime,eval(['vDistance',vStrings{vPlot}]),'Color',vColor);
    end
    
    % prepare the legend
    vNumberOfTracksDisplayed = vNumberOfTracksDisplayed + 1;
    vTrackNames(vNumberOfTracksDisplayed) = cellstr(vTrack2.mName);
end

vStrings = {'Absolute ','X-','Y-','Z-'};
for vPlot = 1:4
    subplot(2,2,vPlot)
    title([vStrings{vPlot}, 'distance from ', vTrack1.mName]);
    legend(vTrackNames);
    xlabel('Time');
    ylabel('Distance [um]');
    hold off
end