%% Written by Q.Peter SU on Jan 28 2016
%% Calculate the STORM SM localization nearestneighbour distance
% based on Rahul's STcoloc.m
% based on Peter's ST09_STcoloc_3D.m
% using filtered precision.mat file output by Yujie's precisionstat.m
%             or    *-Corr.dat file output by DriftCorrect_3D.m
%             or   *Center.mat file output by ST09_2_Convert_Domain_for_Colocalization
%
%% Input
%  1. load *.mat file from precisionstat.m for cluster information
%          *.mat file from Convert_Domain_for_Colocalization.m for cluster after Domain Analysis
%  2. Nearest Neighbour.m
%% Output
%  1. plot nearestneighbour distance histogram with tile
%  2. save two figures




function ST09_3_ST_NearestNeighbour
%% clear
clc
%clear all;
close all
%% 1. Data import
% Red dat file
[fileR pathR] = uigetfile({'*precision*.mat','Precision *.mat file';...
    '*.mat','All mat file';...
    '*filter*.dat','Filtered dots *.dat file';...
    '*.dat','All dat files'},...
    'Choose the Filtered-precision Red Cluster *.mat file.',...
    'Multiselect', 'off');
disp(pathR);
disp('Red DAT File : ');
disp(fileR);
cd(pathR);

% Green dat file
[fileG pathG] = uigetfile({'*precision*.mat','Precision *.mat file';...
    '*.mat','All mat file';...
    '*filter*.dat','Filtered dots *.dat file';...
    '*.dat','All dat files'},...
    'Choose the Filtered-precision Green Cluster *.mat file.',...
    'Multiselect', 'off');
disp('Green DAT File : ');
disp(fileG);
cd(pathG);
red_images = importdata([pathR fileR]);
green_images = importdata([pathG fileG]);
%% 1.1 change matrix
% *.mat file output by Yujie's precisionstat.m
if strcmp(fileR(end-3:end),'.mat');
    % change the red_image to Centers Red
    centersR=[];
    centersR(:,1)=red_images(:,1);   % red cluster center X    in um
    centersR(:,2)=red_images(:,3);   % red cluster center Y    in um
    centersR(:,3)=red_images(:,5);   % red cluster center Z    in um
    % change the green_image to Centers Green
    centersG=[];
    centersG(:,1)=green_images(:,1); % green cluster center X  in um
    centersG(:,2)=green_images(:,3); % green cluster center Y  in um
    centersG(:,3)=green_images(:,5); % green cluster center Z  in um
    
    % *.dat file output by Rahul's Merge_and_plotfile.m and Driftcorrct.m
elseif strcmp(fileR(end-3:end),'.dat');
    % change the red_image to Centers Red
    centersR=[];
    centersR(:,1)=red_images(:,1);   % red cluster center X    in um
    centersR(:,2)=red_images(:,2);   % red cluster center Y    in um
    centersR(:,3)=red_images(:,3);   % red cluster center Z    in um
    % change the green_image to Centers Green
    centersG=[];
    centersG(:,1)=green_images(:,1); % green cluster center X  in um
    centersG(:,2)=green_images(:,2); % green cluster center Y  in um
    centersG(:,3)=green_images(:,3); % green cluster center Z  in um
end

%% 2. input distance between the two color
% if nargin<1
%     distance      =100%input('Enter the distance between the two color in nm. - ');
% end
% distanceRandG=distance/1000;
% if nargin<1
%     dimensionjudge=2;%input('You Want 2D or 3D Colocalization?  Press 2 or 3 - ');
% end

%% 3.1 Nearest neighbour
%% 3.1.1 Calculate the nearestneighbour distance from red to green
%        Finds the nearest neighbour by Euclidean distance to
%        each point in CentersR from CentersG
IDXR2G=nearestneighbour(centersR(:,1:2)',centersG(:,1:2)');
disR2Gxy=sqrt((centersG(IDXR2G,1)-centersR(:,1)).^2 ...
    +(centersG(IDXR2G,2)-centersR(:,2)).^2);
disR2Gz =sqrt((centersG(IDXR2G,3)-centersR(:,3)).^2);
%% 3.1.2 Calculate the nearestneighbour distance from green to red
%        Finds the nearest neighbour by Euclidean distance to
%        each point in CentersG from CentersR
IDXG2R=nearestneighbour(centersG(:,1:2)',centersR(:,1:2)');
disG2Rxy=sqrt((centersR(IDXG2R,1)-centersG(:,1)).^2 ...
    +(centersR(IDXG2R,2)-centersG(:,2)).^2);
disG2Rz =sqrt((centersR(IDXG2R,3)-centersG(:,3)).^2);
%% 4.1 Plot the dots
close all
scrsz = get(0,'ScreenSize');
figure1= figure('OuterPosition',[0.05*scrsz(3) 0.05*scrsz(4) 0.8*scrsz(3) 0.9*scrsz(4)],'Color',...
    [0.313725501298904 0.313725501298904 0.313725501298904]);
subplot(2,2,1)
plot(centersR(:,1),centersR(:,2),'r.');
axis equal
xlabel('X in um');ylabel('Y in um');
title([fileR(1:30)]);

subplot(2,2,2)
plot(centersG(:,1),centersG(:,2),'g.');
axis equal
xlabel('X in um');ylabel('Y in um');
title([fileG(1:30)]);

subplot(2,2,3)
hist(disR2Gxy*1000,20);
xlabel('NearestNeighbour Distance(nm)');ylabel('Occurance');
title({['Red to Green NearestNeighbour Distance'];...
    ['Cluster Centers = ' num2str(size(centersR,1))];...
    ['Mean Distance = ' num2str(mean(disR2Gxy*1000)) 'nm'];...
    ['Median Distance = ' num2str(median(disR2Gxy*1000)) 'nm'];...
    ['std = ' num2str(std(disR2Gxy*1000))]},...
    'FontWeight','bold','FontSize',10,'FontName','Arial',...
    'Color',[1 1 1]);

subplot(2,2,4)
hist(disG2Rxy*1000,20);
xlabel('NearestNeighbour Distance(nm)');ylabel('Occurance');
title({['Green to Red NearestNeighbour Distance'];...
    ['Cluster Centers = ' num2str(size(centersG,1))];...
    ['Mean Distance = ' num2str(mean(disG2Rxy*1000)) 'nm'];...
    ['Median Distance = ' num2str(median(disG2Rxy*1000)) 'nm'];...
    ['std = ' num2str(std(disG2Rxy*1000))]},...
    'FontWeight','bold','FontSize',10,'FontName','Arial',...
    'Color',[1 1 1]);

saveas(gcf, ['0-NearestNeighbour-' fileR(1:25) '.fig'], 'fig');
saveas(gcf, ['0-NearestNeighbour-' fileR(1:25) '.bmp'], 'bmp');



end
