%Hanna, 2005
function [value] = SectorWarpingTool(varargin)
if checkSATDir(pwd) == 0
    fprintf('You are not in a valid SectorAnalysisToolbox Project directory\n');
    return;
end
error(nargchk(0,inf,nargin));
% Open, move, and get the handles from the figure file.
fig = openfig(mfilename, 'reuse');
% Move the gui and then show it, rather than have it jump all over the
% place.
movegui(fig, 'center');
set(fig, 'visible', 'on');
handles = guihandles(fig);
set(fig, 'Color', get(handles.uipanel2, 'BackgroundColor'));


% Set all the callback functions
set(handles.ShowSectorsChk, 'callback', {@doToggleSectors});
set(handles.show_image_chk, 'callback', {@doToggleImage});
set(handles.show_ellipse_chk, 'callback', {@doToggleEllipses});
set(handles.show_sec_num_chk, 'callback', {@doToggleSecNums});
set(handles.calc_mean_shape, 'callback', {@doCalcMeanShape});
set(handles.warp_sectors_btn, 'callback', {@doWarpSectors});
set(handles.ShowEdgeChk, 'callback', {@doToggleEdge});
set(handles.SectorFilesPopup, 'callback', {@doSkipFile});
set(handles.nextBtn, 'callback', {@doNext});
set(handles.prevBtn, 'callback', {@doPrev});
set(handles.zoom_btn, 'callback', {@doZoom});
set(handles.pan_btn, 'callback', {@doPan});
set(handles.reorderEdge_btn, 'callback', {@doClickControlPoints});
set(handles.original_pos_radio, 'callback', {@switch_sectors});
set(handles.aligned_pos_radio, 'callback', {@switch_sectors});
set(handles.warped_pos_radio, 'callback', {@switch_sectors});
set(handles.original_edge_radio, 'callback', {@switch_edges});
set(handles.aligned_edge_radio, 'callback', {@switch_edges});
set(handles.mean_edge_radio, 'callback', {@switch_edges});
set(handles.flip_lr_btn, 'callback', {@doFlipLR});
set(handles.reverse_ordering_btn, 'callback', {@doReverseOrdering});
set(handles.calc_ellipse_info_btn, 'callback', {@doCalcRegionProps});
set(handles.imageAxis, 'ButtonDownFcn', {@clickedInWindow});
set(handles.meanShapeScaleNorm, 'callback', {@doToggleMeanShapeCheckboxs});
set(handles.scaleOpt, 'callback', {@doToggleScaleopt});


%set(handles.move_edge_pt_btn, 'callback', {@doMoveEdgePoints});


% Initialize the application data structure
ad.figMain = fig;

ad.StagedDirString = 'Staged';
ad.OriginalsDirString = 'Originals';
ad.ProcessedDirString = 'Processed';
ad.DataDirString = [ad.StagedDirString, filesep, 'Data'];
if ~exist(ad.DataDirString, 'dir'); mkdir(ad.DataDirString); end;
ad.projectpath = pwd;
ad.sector_files = FindFilenames('_sectors.tif', ad.StagedDirString);
if isempty(ad.sector_files)
    delete(fig);
    return;
end
if ~exist([ad.DataDirString, filesep, 'mean_shape.mat'], 'file')
    file = ad.sector_files{1};
    [path, file, ext, vers] = fileparts(file);
    ad.mean_shape = load([path, filesep, 'edge_xy.mat']);
    ad.mean_shape = ad.mean_shape.E;
    mean_shape = ad.mean_shape;
    save([ad.DataDirString, filesep, 'mean_shape.mat'], 'mean_shape');
    clear mean_shape;
else
    ad.mean_shape = load([ad.DataDirString, filesep, 'mean_shape.mat']);
    ad.mean_shape = ad.mean_shape.mean_shape;
end

ad.value = 0;
ad.handles = handles;
ad.sectors_handle =  [];
ad.ellipse_handle =  [];
ad.mean_shape_handle = [];
ad.sector_im = [];
ad.ellipse_axis_handle = [];
ad.aligned_edge_handle = [];
ad.aligned_sectors_handle = [];
ad.warped_sectors_handle = [];
ad.control_marker_handle = [];
ad.edge_handle  = [];
ad.text_edge_handle = [];
ad.aligned_text_edge_handle = [];
ad.mean_text_edge_handle = [];
ad.image_handle = [];
ad.marker_indices = 1;
ad.edge_pt_to_move = [];
setappdata(0,'SectorWarpingToolData',ad);
init_gui(ad, handles); % there is an error in here;

set(ad.figMain,'visible','on')
% try
%     uiwait(fig);
% catch
%     if ishandle(fig)
%         delete(fig)
%     end
% end
value = 0;
return;

function [] = init_gui(ad, handles)
if isempty(ad.sector_files)
    str = 'No sector files found';
    set(handles.SectorFilesPopup, 'String', str);
    return;
end
str = ad.sector_files;
set(handles.SectorFilesPopup, 'String', str);
name = str{1};
loadImage(name);
loadSectors(name);
loadEdges(name);
displayImage();
displayEdges();
displaySectors();
displayEllipses();

return;
%%%
%
%%%
function ad = init_handles(ad)
ad.sectors_handle =  [];
ad.ellipse_axis_handle = [];
ad.ellipse_handle =  [];
ad.mean_shape_handle = [];
ad.aligned_edge_handle = [];
ad.aligned_sectors_handle = [];
ad.warped_sectors_handle = [];
ad.control_marker_handle = [];
ad.edge_handle  = [];
ad.image_handle = [];
ad.edge_pt_to_move = [];
return
%%%
%
%%%
function doSkipFile(sector_popup, evd)
ad =getappdata(0, 'SectorWarpingToolData');
set(ad.figMain, 'Pointer', 'watch');
if ishandle(ad.ellipse_handle)
delete(ad.ellipse_handle);
end
cla(ad.handles.imageAxis);
ad = init_handles(ad);
if isempty(ad.sector_files); return; end;
zoom off; pan off;
set(ad.handles.zoom_btn, 'value', 0);
set(ad.handles.pan_btn, 'value', 0);
value = get(sector_popup, 'value');
fprintf('Skipped to file %s\n', ad.sector_files{value});
setappdata(0,'SectorWarpingToolData',ad);
loadImage(ad.sector_files{value});
loadSectors(ad.sector_files{value});
loadEdges(ad.sector_files{value});
displayImage();
displayEdges();
displaySectors();
displayEllipses();
set(ad.figMain, 'Pointer', 'arrow');

return;

function doFlipLR(flip_lr_btn, evd)
ad =getappdata(0, 'SectorWarpingToolData');
set(ad.figMain, 'Pointer', 'watch');
s = get(ad.handles.SectorFilesPopup, 'String');
v = get(ad.handles.SectorFilesPopup, 'Value');
file = s{v};
[path, file, ext, vers] = fileparts(file);
file = file(1:end-8);
organ_filename = [path, filesep, file, '_organ.tif'];
sectors_filename = [path, filesep, file, '_sectors.tif'];
image_filename = [path, filesep, file, '.JPG'];
if exist(organ_filename, 'file')
    I = imread(organ_filename);
    I = fliplr(I);
    imwrite(I, organ_filename, 'TIFF');
    E = binim2SortedEdge(I);
    save([path, filesep, 'edge_xy.mat'], 'E');

end
if exist(sectors_filename, 'file')
    I = imread(sectors_filename);
    I = fliplr(I);
    imwrite(I, sectors_filename, 'TIFF');
    [x, y] = find(I >0);
    sectors = [x(:), y(:)];
    save([path, filesep, 'sectors_xy.mat'], 'sectors');
    sector_info = sector_im2struct(I);
    save([path, filesep, 'sector_info.mat'], 'sector_info');
end
if exist(image_filename, 'file')
    I = imread(image_filename);
    I = flipRGB(I);
    imwrite(I, image_filename, 'JPEG');
end
doSkipFile(ad.handles.SectorFilesPopup);
set(ad.figMain, 'Pointer', 'arrow');
setappdata(0, 'SectorWarpingToolData', ad);
%%%%
%
%
%%%%
% function doMoveEdgePoints(move_edge_pt_btn, evd)
% ad =getappdata(0, 'SectorWarpingToolData');
% switch (get(move_edge_pt_btn, 'Value'))
%     case 0
%         set(ad.figMain, 'Pointer', 'Arrow');
%         set(ad.figMain, 'WindowButtonDownFcn', '', 'WindowButtonMotionFcn' , '', 'WindowButtonUpFcn', '', 'BackingStore','on');
%         set(ad.handles.imageAxis, 'ButtonDownFcn', {@doClickControlPoints});
%         filename = get(ad.handles.SectorFilesPopup, 'string');
%         filename = filename{get(ad.handles.SectorFilesPopup, 'value')};
%         [pathname, filename, ext, vers] = fileparts(filename);
%         E = ad.edge;
%         save([pathname, filesep, 'edge_xy.mat'], 'E');
%     case 1
%         set(ad.handles.original_edge_radio, 'value', 1);
%         setappdata(0, 'SectorWarpingToolData', ad);
%         switch_edges;
%         set(ad.figMain, 'Pointer', 'crosshair');
%         set(ad.handles.imageAxis, 'ButtonDownFcn', '');
%         set(ad.figMain, 'BackingStore', 'off', 'WindowButtonDownFcn',  {@select_edge_pt});
% end
% setappdata(0, 'SectorWarpingToolData', ad);
% return;
%%%%
%
%
%
%%%%
% function select_edge_pt(fig, evd)
% ad =getappdata(0, 'SectorWarpingToolData');
% cp = get(ad.handles.imageAxis, 'CurrentPoint');
% cp = cp(1,1:2);
% [val, ind] = min(sqrt(sum((fliplr(ad.edge) - ones(size(ad.edge,1),1)*cp).^2,2)));
% ad.edge_pt_to_move = ind(1);
% setappdata(0, 'SectorWarpingToolData', ad);
% set(fig, 'WindowButtonMotionFcn', {@move_edge_pt}, 'WindowButtonUpFcn',  {@drop_edge_pt});
% return;
%%%%
%
%
%
%%%%
% function move_edge_pt(fig, evd)
% ad =getappdata(0, 'SectorWarpingToolData');
% if isempty(ad.edge_pt_to_move)
%     return;
% end
% cp = get(ad.handles.imageAxis, 'CurrentPoint');
% cp = cp(1,1:2);
% ad.edge(ad.edge_pt_to_move(1),:) = fliplr(cp);
% if length(ad.edge_handle)== 0
%     return
% end
% ad.edge_handle = ad.edge_handle(1);
% xdata = get(ad.edge_handle, 'XData');
% ydata = get(ad.edge_handle, 'YData');
% if iscell(xdata)
%     xdata = xdata{1};
%     ydata = ydata{1};
% end
% xdata(ad.edge_pt_to_move(1)) = cp(1);
% ydata(ad.edge_pt_to_move(1)) = cp(2);
% set(ad.edge_handle, 'XData', xdata);
% set(ad.edge_handle, 'YData', ydata);
% setappdata(0, 'SectorWarpingToolData', ad);
% return;
%%%%
%
%
%
%%%%
% function drop_edge_pt(fig, evd)
% ad =getappdata(0, 'SectorWarpingToolData');
% set(fig, 'WindowButtonMotionFcn' ,'','WindowButtonUpFcn','');
% ad.ctrl_pt_to_move = [];
% setappdata(0, 'SectorWarpingToolData', ad);
% return;
%%%
%
%%%
function doNext(next_btn, evd)

ad =getappdata(0, 'SectorWarpingToolData');
if isempty(ad.sector_files); return; end;
zoom off; pan off;
set(ad.handles.zoom_btn, 'value', 0);
set(ad.handles.pan_btn, 'value', 0);
value = get(ad.handles.SectorFilesPopup, 'value');
if value < length(ad.sector_files)
    value = value + 1;
    set(ad.handles.SectorFilesPopup, 'value', value);
    doSkipFile(ad.handles.SectorFilesPopup);
end
return;


function doPrev(prev_btn, evd)
ad =getappdata(0, 'SectorWarpingToolData');
if isempty(ad.sector_files); return; end;
value = get(ad.handles.SectorFilesPopup, 'value');
if value > 1
    value = value - 1;
    set(ad.handles.SectorFilesPopup, 'value', value);
    doSkipFile(ad.handles.SectorFilesPopup);
end
return;
%%%
%
%%%
% function doToggleMeanShape(mean_shape_chk, evd)
% ad =getappdata(0, 'SectorWarpingToolData');
% if get(mean_shape_chk, 'value')
%     set(ad.mean_shape_handle, 'Visible', 'on');
% else
%     set(ad.mean_shape_handle, 'Visible', 'off');
% end
% setappdata(0,'SectorWarpingToolData',ad);
%%
%
function switch_sectors(toggle_sector_chk, evd)
ad =getappdata(0, 'SectorWarpingToolData');
if get(ad.handles.original_pos_radio, 'value')
    set(ad.sectors_handle, 'Visible', 'on');
    set(ad.warped_sectors_handle, 'Visible', 'off');
    set(ad.aligned_sectors_handle, 'Visible', 'off');

elseif get(ad.handles.aligned_pos_radio, 'value')
    set(ad.aligned_sectors_handle, 'Visible', 'on');
    set(ad.sectors_handle, 'Visible', 'off');
    set(ad.warped_sectors_handle, 'Visible', 'off');
else
    set(ad.aligned_sectors_handle, 'Visible', 'off');
    set(ad.sectors_handle, 'Visible', 'off');
    set(ad.warped_sectors_handle, 'Visible', 'on');
end
setappdata(0,'SectorWarpingToolData',ad);
displayEllipses;
%%
%
function switch_edges(toggle_sector_chk, evd)
ad =getappdata(0, 'SectorWarpingToolData');
if get(ad.handles.original_edge_radio, 'value')
    set(ad.edge_handle, 'Visible', 'on');
    set(ad.aligned_edge_handle, 'Visible', 'off');
    set(ad.mean_shape_handle, 'Visible', 'off');
    set(ad.aligned_text_edge_handle, 'Visible', 'off');
    set(ad.mean_text_edge_handle, 'Visible', 'off');
    set(ad.text_edge_handle, 'Visible', 'on');
elseif get(ad.handles.aligned_edge_radio, 'value')
    set(ad.edge_handle, 'Visible', 'off');
    set(ad.aligned_edge_handle, 'Visible', 'on');
    set(ad.mean_shape_handle, 'Visible', 'off');
    set(ad.aligned_text_edge_handle, 'Visible', 'on');
    set(ad.mean_text_edge_handle, 'Visible', 'off');
    set(ad.text_edge_handle, 'Visible', 'off');
else
    set(ad.edge_handle, 'Visible', 'off');
    set(ad.aligned_edge_handle, 'Visible', 'off');
    set(ad.mean_shape_handle, 'Visible', 'on');
    set(ad.aligned_text_edge_handle, 'Visible', 'off');
    set(ad.mean_text_edge_handle, 'Visible', 'on');
    set(ad.text_edge_handle, 'Visible', 'off');
end

function doToggleEdge(toggle_edge_chk, evd)
ad =getappdata(0, 'SectorWarpingToolData');
if get(toggle_edge_chk, 'value')
    set(ad.handles.original_edge_radio, 'Enable', 'on');
    set(ad.handles.aligned_edge_radio, 'Enable', 'on');
    set(ad.handles.mean_edge_radio, 'Enable', 'on');
    switch_edges;
else
    set(ad.handles.original_edge_radio, 'Enable', 'off');
    set(ad.handles.aligned_edge_radio, 'Enable', 'off');
    set(ad.handles.mean_edge_radio, 'Enable', 'off');
    set(ad.edge_handle, 'Visible', 'off');
    set(ad.aligned_edge_handle, 'Visible', 'off');
    set(ad.mean_shape_handle, 'Visible', 'off');

end
setappdata(0,'SectorWarpingToolData',ad);


function doToggleSectors(toggle_sector_chk, evd)
ad =getappdata(0, 'SectorWarpingToolData');
if get(toggle_sector_chk, 'value')
    set(ad.handles.original_pos_radio, 'Enable', 'on');
    set(ad.handles.aligned_pos_radio, 'Enable', 'on');
    set(ad.handles.warped_pos_radio, 'Enable', 'on');
    switch_sectors;
    ad =getappdata(0, 'SectorWarpingToolData');

else
    set(ad.sectors_handle, 'Visible', 'off');
    set(ad.aligned_sectors_handle, 'Visible', 'off');
    set(ad.warped_sectors_handle, 'Visible', 'off');
    set(ad.handles.original_pos_radio, 'Enable', 'off');
    set(ad.handles.aligned_pos_radio, 'Enable', 'off');
    set(ad.handles.warped_pos_radio, 'Enable', 'off');
end
setappdata(0,'SectorWarpingToolData',ad);
%%%
%
%%%
function doZoom(zoom_btn, evd)
ad =getappdata(0, 'SectorWarpingToolData');
zoom off;
pan off;
if get(ad.handles.zoom_btn, 'value')
    zoom on;
end
%%%
%
%%%
% function doToggleAll()
% ad =getappdata(0, 'SectorWarpingToolData');
% ishandle(ad.ellipse_handle)
% doToggleSectors(ad.handles.ShowSectorsChk);
% doToggleImage(ad.handles.show_image_chk);
% doToggleEllipses(ad.handles.show_ellipse_chk);
% doToggleEdge(ad.handles.ShowEdgeChk);
% return;

%%%
%
%%%
% function doAlignMeanShape(align_mean_shape_btn, evd)
% ad =getappdata(0, 'SectorWarpingToolData');
% ms = ad.mean_shape;
% set(ad.figMain, 'Visible', 'off');
% [mean_shape, val] = AlignShape('Shape', ms);
% ad.mean_shape = mean_shape;
% save([ad.DataDirString, filesep, 'mean_shape.mat'], 'mean_shape');
% set(ad.figMain, 'Visible', 'on');
% setappdata(0, 'SectorWarpingToolData', ad);
% updateGraphicalData;
%%%
%
%%%
% function dropshape()
% ad =getappdata(0, 'SectorWarpingToolData');
% set(ad.figMain, 'WindowButtonDownFcn', '', 'WindowButtonDownFcn', '');
% delete(ad.ms_handle);
% % Put the settings back as they were.
% set(ad.handles.ShowSectorsChk, 'Value', show_sec_val);
% set(ad.handles.ShowEdgeChk, 'Value', show_edge_val);
% set(ad.handles.show_image_chk, 'Value', show_image_val);
% set(ad.handles.show_ellipse_chk, 'Value', show_ellip_val);
% doToggleAll
% setappdata(0, 'SectorWarpingToolData', ad);
%%%
%
%%%
function doPan(pan_btn, evd)
ad =getappdata(0, 'SectorWarpingToolData');
zoom off;
pan off;
if get(ad.handles.pan_btn, 'value')
    pan on;
end
%%%
%
%%%
function doClickControlPoints(image_axis, evd)

ad =getappdata(0, 'SectorWarpingToolData');
if strcmp(get(ad.figMain, 'SelectionType'), 'alt')
    doNext;
else
    %c = get(image_axis, 'CurrentPoint');    
    %x = c(1,1); y = c(1, 2);
    
    [x,y] = ginput(1);
    d = sqrt((ad.edge(:,2) - x).^2 + (ad.edge(:,1) - y).^2);
    [v, i] = min(d);
    i = i(1);
    e_temp = ad.edge;
    ad.edge = [e_temp(i:end,:); e_temp(1:i-1,:)];
    E = ad.edge;
    value = get(ad.handles.SectorFilesPopup, 'value');
    disp(ad.sector_files{value});
    [path, name, ext, vers] = fileparts(ad.sector_files{value});
    save([path, filesep, 'edge_xy.mat'], 'E');
    setappdata(0,'SectorWarpingToolData',ad);
    displayEdges();
    displaySectors();
end
%%%
%
%%%
function doReverseOrdering(doreverse_order_btn, evd)
ad =getappdata(0, 'SectorWarpingToolData');
e_temp = ad.edge;
ad.edge = flipud(e_temp);
E = ad.edge;
value = get(ad.handles.SectorFilesPopup, 'value');
disp(ad.sector_files{value});
[path, name, ext, vers] = fileparts(ad.sector_files{value});
save([path, filesep, 'edge_xy.mat'], 'E');
setappdata(0,'SectorWarpingToolData',ad);
displayEdges();
displaySectors();
%%%
%
%%%



function doCalcRegionProps(calc_region_props_btn, evd) %#ok<INUSD,INUSD>
ad =getappdata(0, 'SectorWarpingToolData');
StagedDirString = ad.StagedDirString;
FILES = findfilesext('jpg', StagedDirString);
if isempty(FILES); return; end;
[s,v] = listdlg('Name', 'Sectors To Warp', 'PromptString','Select the files you wish to calc. ellipse information for:', 'SelectionMode','multiple', 'ListString',FILES, 'ListSize', [600, 300]);
if isempty(s)
    fprintf('User did not select any files\n');
    return;
end
files = FILES(s);
ellipfitalg=questdlg('Choose Fitting Algorithm?', 'Ellipse Fitting Algorithm', 'AG', 'AIH', 'EEK','AIH');

h = waitbar(0, 'Calculating ellipse information...');
for i=1:length(files)
    file = files{i};
    fprintf('Calculating ellipse in image %s\n', file);
    [path, file, ext, vers] = fileparts(file);
    
    % Calculate the ellipse information for the warped sectors
    cellplanes = [path, filesep, [file, '_sectors.tif']];
    warped_filename = [path, filesep, 'warped_sectors_xy.mat'];
    if exist(warped_filename, 'file')
        warped_sectors = load(warped_filename);
        warped_sectors = warped_sectors.warped_sectors;
        warped_sector_info = sector_im2struct(warped_sectors, ellipfitalg, cellplanes);
        save([path, filesep, 'warped_sector_info.mat'], 'warped_sector_info');
        origname = [path, filesep, 'warped_sector_info_orig.mat'];
        if exist(origname, 'file')
            delete(origname);
        end
    end
    waitbar(i/length(files), h);
end
close(h);
fprintf('done');
updateGraphicalData;
%%%
%
%%%
function doWarpSectors(calc_mean_shape_btn, evd) %#ok<INUSD,INUSD>
ad =getappdata(0, 'SectorWarpingToolData');
StagedDirString = ad.StagedDirString;
FILES = findfilesext('jpg', StagedDirString);
if isempty(FILES); return; end;
[s,v] = listdlg('Name', 'Sectors To Warp', 'PromptString','Select the files you wish to warp to the mean shape:', 'SelectionMode','multiple', 'ListString',FILES, 'ListSize', [500, 200]);
if v == 0
    return;
end
files = FILES(s);
mean_shape = load([ad.DataDirString, filesep, 'mean_shape.mat']);
mean_shape = mean_shape.mean_shape;


warpList{1} = 'Petiole and Lamina (Sam)';
warpList{2} = 'Concave Shape (Erika)';
warpList{3} = 'Global Warp (Nobody)';
warpList{4} = 'Petals - Concave Shape (Sue)';
[warpSelection,v] = listdlg('PromptString','Select a warp:',...
                'SelectionMode','single',...
                'ListString',warpList);
 

h = waitbar(0, 'Warping sectors please wait...');
for i=1:length(files)
    file = files{i}
    fprintf('Warping sectors in image %s\n', file);
    [path, file, ext, vers] = fileparts(file);
    if ~exist([path, filesep, 'edge_aligned_xy.mat'], 'file')
        fprintf('The edge for %s has not been aligned', [path, filesep, name]);
    else
        edge = load([path, filesep, 'edge_xy.mat']);
        edge = edge.E;
        aligned_edge = load([path, filesep, 'edge_aligned_xy.mat']);
        aligned_edge = aligned_edge.edge_aligned_xy;
        sectors = imread([path, filesep, [file, '_sectors.tif']]);
        sz = size(sectors);
        boundaries = bwboundaries(sectors);
        if(warpSelection == 1)%sam warp
            aligned_sectors = warpSectorsViaTriangles(aligned_edge, edge, boundaries);
            warped_sectors = warpSectorsViaTriangles(mean_shape, aligned_edge, aligned_sectors);
        elseif(warpSelection == 2)%erika warp
            aligned_sectors = warpSectorsViaTrianglesConcave(aligned_edge, edge, boundaries);
            warped_sectors = warpSectorsViaTrianglesConcave(mean_shape, aligned_edge, aligned_sectors);
        elseif(warpSelection == 3)%old warp
            aligned_sectors = warp_sectors(aligned_edge, edge, boundaries);
            warped_sectors = warp_sectors(mean_shape, aligned_edge,aligned_sectors);
        elseif(warpSelection == 4)%erika warp
            aligned_sectors = warpSectorsViaTrianglesConcave(aligned_edge, edge, boundaries);
            warped_sectors = warpSectorsViaTrianglesConcave(mean_shape, aligned_edge, aligned_sectors);
        else
            error('No warpSelection')            
        end
          
        save([path, filesep, 'aligned_sectors_xy.mat'], 'aligned_sectors');
        save([path, filesep, 'warped_sectors_xy.mat'], 'warped_sectors');
        
        %%% code for figure to show the difference between global and pwl
        %%% warp
        figure;
        %plot(edge(:,1),edge(:,2),'rx');
        hold on
        plot(aligned_edge(:,1),aligned_edge(:,2),'bo')
        for j = 1:numel(aligned_sectors)
            as = aligned_sectors{j};
            plot(as(:,1),as(:,2),'bx')
        end
        plot(mean_shape(:,1),mean_shape(:,2),'go');
        for j = 1:size(mean_shape,1)
         %   text(mean_shape(j,1),mean_shape(j,2),num2str(j))
        end
        
        
        for j = 1:numel(warped_sectors)
            ws = warped_sectors{j};
            plot(ws(:,1),ws(:,2),'gx')
            
        end
        title(path);
        
        paul_warped_sectors = warp_sectors(mean_shape, edge, boundaries);
        for j = 1:numel(paul_warped_sectors)
            pws = paul_warped_sectors{j};
            plot(pws(:,1),pws(:,2),'yx')
        end
        %%%
        
        
        
        
    end
    waitbar(i/length(files), h);
end
close(h);
updateGraphicalData;

%%%
% function warped_sectors = warp_sectors(warped_shape, orig_shape, sectors)
% Tinv = cp2tform(orig_shape, warped_shape, 'linear conformal');
% warped_sectors = sectors;
% for b=1:length(sectors)
%     s = sectors{b};
%     x =s(:,1); y = s(:,2);
%     [xm, ym] = tformfwd(Tinv, x, y);  % Results should equal [x, y]
%     warped_sectors{b} = [xm(:), ym(:)];
% end
% return;
%%%
%
%%%

function warped_sectors = warp_sectors(warped_shape, orig_shape, sectors)
Tinv = cp2tform(warped_shape, orig_shape, 'polynomial', 2);
warped_sectors = sectors;
for b=1:length(sectors)
    s = sectors{b};
    x =s(:,1); y = s(:,2);
    warped_sectors{b} = round([ones(length(x),1)  x  y  x.*y  x.^2  y.^2]*Tinv.tdata);
end
return;

% function warped_sectors = warp_sectors(warped_shape, orig_shape, sectors)
% Tinv = cp2tform(warped_shape, orig_shape, 'polynomial', 3);
% warped_sectors = sectors;
% for b=1:length(sectors)
%     s = sectors{b};
%     x =s(:,1); y = s(:,2);
%     warped_sectors{b} = round([ones(length(x),1) x  y  x.*y  x.^2  y.^2  y.*x.^2  x.*y.^2  x.^3  y.^3]*Tinv.tdata);
% end
% return;



%%%%
%
%
%%%%
function updateGraphicalData()
ad =getappdata(0, 'SectorWarpingToolData');
loadImage(ad.sector_files{get(ad.handles.SectorFilesPopup, 'value')});
loadSectors(ad.sector_files{get(ad.handles.SectorFilesPopup, 'value')});
loadEdges(ad.sector_files{get(ad.handles.SectorFilesPopup, 'value')});
displayEdges();
displaySectors();
%%%%
%
%
%%%%


function doToggleMeanShapeCheckboxs(meanShapeScaleNorm,evd)
ad =getappdata(0, 'SectorWarpingToolData');
if(get(ad.handles.meanShapeScaleNorm,'Value') == 1)
    set(ad.handles.scaleOpt,'Value',1)
end

return

function doToggleScaleopt(scaleOpt,evd)
ad =getappdata(0, 'SectorWarpingToolData');
if(get(ad.handles.meanShapeScaleNorm,'Value') == 1)
    set(ad.handles.scaleOpt,'Value',1)
end
return

function doCalcMeanShape(calc_mean_shape_btn, evd) %#ok<INUSD,INUSD>
ad =getappdata(0, 'SectorWarpingToolData');
files = ad.sector_files;
x = [];
h = waitbar(0, 'Calculating the mean organ shape please wait...');
for f = 1:length(files)
    waitbar(f/length(files), h);
    file = files{f};
    [path, file, ext, vers] = fileparts(file);
    if ~exist([path, filesep, 'edge_xy.mat'], 'file')
        fprintf('There is no edge information for %s\n', path);
        return;
    else
        edge = load([path, filesep, 'edge_xy.mat']);
        edge = edge.E';
        x = cat(2, x, edge(:));
    end
end
close(h);
opts.scaling = get(ad.handles.scaleOpt,'Value');
opts.translation = 1;
opts.rotation = 1;
opts.rescaling = get(ad.handles.meanShapeScaleNorm,'Value');


fid = fopen([ad.DataDirString, filesep,'ProcrustesOptions.txt'],'w');
fprintf(fid,'scaling:%6.0f\n',opts.scaling);
fprintf(fid,'translation:%6.0f\n',opts.translation);
fprintf(fid,'rotation:%6.0f\n',opts.rotation);
fprintf(fid,'re-scaling:%6.0f\n',opts.rescaling);
fclose(fid);


val = get(ad.handles.meanShapeScaleNorm,'Value');
if(val)
    [aligned_pts, template] = pcalib_GPAMeanScale('data', x, 'opts', opts);
else
    [aligned_pts, template] = pcalib_GPA('data', x, 'opts', opts);
end

%an intereting line of commented code below
%[aligned_pts, template] = pcalib_GPA('data', x, 'opts', opts, 'template', mean(x,2));

mean_shape = template';
ad.mean_shape = mean_shape;
save([ad.DataDirString, filesep, 'mean_shape.mat'], 'mean_shape');
for f = 1:length(files)
    file = files{f};
    [path, file, ext, vers] = fileparts(file);
    if ~exist([path, filesep, 'edge_xy.mat'], 'file')
    else
        apts = aligned_pts(:,f);
        edge_aligned_xy = (reshape(apts, 2, length(apts)/2))';
        save([path, filesep, 'edge_aligned_xy.mat'], 'edge_aligned_xy');
    end
end
setappdata(0,'SectorWarpingToolData',ad);
updateGraphicalData;
disp('Calculating the mean shape');
%%%%
%
%
%%%%
function doToggleEllipses(show_ellipse_chk, evd)
ad =getappdata(0, 'SectorWarpingToolData');
if ishandle(ad.ellipse_handle)
    if get(show_ellipse_chk, 'value')
        set(ad.ellipse_handle, 'Visible', 'on');
        set(ad.ellipse_axis_handle, 'Visible', 'on');
    else
        set(ad.ellipse_handle, 'Visible', 'off');
        set(ad.ellipse_axis_handle, 'Visible', 'off');
    end
end
setappdata(0,'SectorWarpingToolData',ad);
%%%%
%
%
%%%%
function doToggleSecNums(show_sec_num_chk, evd)
figure;
return
%%%%
%
%
%%%%
function doToggleImage(show_image_btn, evd)
ad =getappdata(0, 'SectorWarpingToolData');
if get(show_image_btn, 'value')
    if isempty(ad.image_handle)
        value = get(ad.handles.SectorFilesPopup, 'value');
        loadImage(ad.sector_files{value});
        ad =getappdata(0, 'SectorWarpingToolData');
        displayImage;
        displayEdges();
        displaySectors();
        displayEllipses();
    end
    set(ad.image_handle, 'Visible', 'on');
else
    set(ad.image_handle, 'Visible', 'off');
end
setappdata(0,'SectorWarpingToolData',ad);
%%%%
%
%
%%%%
function loadSectors(name)
[path, name, ext, vers] = fileparts(name);
ad.sector.x = [];
ad.sector.y = [];
ad =getappdata(0, 'SectorWarpingToolData');
if ~exist([path, filesep, 'sectors_xy.mat'], 'file')
    fprintf('There is no sector information for %s\n', path);
    return;
end
sectors = load([path, filesep, 'sectors_xy.mat']);
ad.sector.x = sectors.sectors(:,1);
ad.sector.y = sectors.sectors(:,2);
if ~exist([path, filesep, 'aligned_sectors_xy.mat'], 'file')
    ad.aligned_sector = ad.sector;
else
    aligned_sectors = load([path, filesep, 'aligned_sectors_xy.mat']);
    aligned_sectors = aligned_sectors.aligned_sectors;

    x =[]; y = [];
        if iscell(aligned_sectors)
    for i=1:length(aligned_sectors); x = [x; aligned_sectors{i}(:,1)]; y = [y; aligned_sectors{i}(:,2)]; end
        end
    ad.aligned_sector.x = x;
    ad.aligned_sector.y = y;

    
end
if ~exist([path, filesep, 'warped_sectors_xy.mat'], 'file')
    ad.warped_sector = ad.sector;
else
    warped_sectors = load([path, filesep, 'warped_sectors_xy.mat']);
    warped_sectors = warped_sectors.warped_sectors;

     x =[]; y = [];
         if iscell(warped_sectors)
    for i=1:length(warped_sectors); x = [x; warped_sectors{i}(:,1)]; y = [y; warped_sectors{i}(:,2)]; end
         end
    ad.warped_sector.x = x;
    ad.warped_sector.y = y;
end
setappdata(0,'SectorWarpingToolData',ad);
return;
%%%%
%
%
%%%%
function loadEdges(name)
[path, name, ext, vers] = fileparts(name);
ad =getappdata(0, 'SectorWarpingToolData');
if ~exist([path, filesep, 'edge_xy.mat'], 'file')
    fprintf('There is no edge information for %s\n', path);
    return;
end
ad.edge = load([path, filesep, 'edge_xy.mat']);
ad.edge = ad.edge.E;
if ~exist([path, filesep, 'edge_aligned_xy.mat'], 'file')
    ad.aligned_edge = ad.edge;
else
    ad.aligned_edge = load([path, filesep, 'edge_aligned_xy.mat']);
    ad.aligned_edge = ad.aligned_edge.edge_aligned_xy;
end
setappdata(0,'SectorWarpingToolData',ad);
return;
%%%
%
%%%
function [sector_im, sector_info] = calcOriginalTensors(pathname, filename)
sector_im = imread([pathname, filesep, filename, '.tif']);
[x,y] = find(sector_im > 0);
sectors = [x(:),y(:)];
save([pathname, filesep, 'sectors_xy.mat'],'sectors');
boundaries = bwboundaries(sector_im);
sector_info = sector_im2struct(boundaries, 'AIH');
save([pathname, filesep, 'sector_info.mat'], 'sector_info');
return
%%%%
%
%
%%%%
function loadImage(name)
ad =getappdata(0, 'SectorWarpingToolData');
[path, name, ext, vers] = fileparts(name);
if get(ad.handles.show_image_chk, 'value')
    ad.currentImage = imread([path, filesep, name(1:end-8), '.JPG']);
else
    ad.currentImage = [];
end
if ~exist([path, filesep, 'sector_info.mat'], 'file')
    [ad.sector_im, ad.sector_info] = calcOriginalTensors(path, name);
else
    ad.sector_info = load([path, filesep, 'sector_info.mat']);
    ad.sector_info = ad.sector_info.sector_info;
    if ~isfield(ad.sector_info(1), 'growthTensor')
        [ad.sector_im, ad.sector_info] = calcOriginalTensors(path, name);
    end
end
if ~exist([path, filesep, 'aligned_sector_info.mat'], 'file')
    ad.aligned_sector_info = ad.sector_info;
else
    ad.aligned_sector_info = load([path, filesep, 'aligned_sector_info.mat']);
    ad.aligned_sector_info = ad.aligned_sector_info.aligned_sector_info;
end
if ~exist([path, filesep, 'warped_sector_info.mat'], 'file')
    ad.warped_sector_info = ad.sector_info;
else
    ad.warped_sector_info = load([path, filesep, 'warped_sector_info.mat']);
    ad.warped_sector_info = ad.warped_sector_info.warped_sector_info;
end
setappdata(0,'SectorWarpingToolData',ad);
return;
%%%%
%
%
%%%%
function displayImage()
ad =getappdata(0, 'SectorWarpingToolData');
if ishandle(ad.image_handle)
    delete(ad.image_handle);
end
if ~isempty(ad.currentImage) && (get(ad.handles.show_image_chk, 'value'))
    %value = get(ad.handles.SectorFilesPopup, 'value');
    %loadImage(ad.sector_files{value});

    %ad =getappdata(0, 'SectorWarpingToolData');

    ad.image_handle = imagesc(ad.currentImage, 'Parent', ad.handles.imageAxis, 'HitTest', 'off');
    hold(ad.handles.imageAxis, 'on');
    axis(ad.handles.imageAxis, 'image', 'ij');
    setappdata(0,'SectorWarpingToolData',ad);
    doToggleImage(ad.handles.show_image_chk);
end
%%%%
%
%
%%%%
function displayEllipses()
ad =getappdata(0, 'SectorWarpingToolData');
if ishandle(ad.ellipse_handle)
    delete(ad.ellipse_handle);
end
if ishandle(ad.ellipse_axis_handle)
    delete(ad.ellipse_axis_handle);
end
ad.ellipse_axis_handle = [];
ad.ellipse_handle = [];
if get(ad.handles.original_pos_radio, 'value')
    sector_info = ad.sector_info;
elseif get(ad.handles.aligned_pos_radio, 'value')
    sector_info = ad.aligned_sector_info;
else
    sector_info = ad.warped_sector_info;
end
try
for i=1:length(sector_info)
    r = sector_info(i);
    ad.ellipse_handle = [ad.ellipse_handle; plot_ellipse(ad.handles.imageAxis, r)];
end
catch
    fprintf('No growth tensor information available. Just showing sector/division plane outlines\n');
    if isempty(axis_handle)
    axis_handle = gca;
    end
end
   % [path, file, ext, vers] = fileparts(file);
  %  ad.mean_shape = load([path, filesep, 'edge_xy.mat']);
%%% load display 'aligned_sectors_xy.mat'
%'warped_sectors_xy.mat'
%    hold(axis_handle, 'on');

 
    %plot(axis_handle, Y(1, :) + offset(1), Y(2,:) + offset(2), 'Color', col, 'HitTest', 'off', 'LineWidth', linewidth);
  %  cross_ph = gtlib_plotGrowthTensorCross('growth_tensor', GT, 'offset', offset, 'colour', col, 'parent', axis_handle, 'scale', scale, 'linewidth', linewidth);


set(ad.ellipse_handle, 'HitTest', 'off');
%set(ad.handles.show_ellipse_chk, 'value', 1);
setappdata(0,'SectorWarpingToolData',ad);
doToggleEllipses(ad.handles.show_ellipse_chk);
setappdata(0,'SectorWarpingToolData',ad);
%%%%
%
%
%%%%
function displayEdges()
ad =getappdata(0, 'SectorWarpingToolData');
if ishandle(ad.edge_handle); delete(ad.edge_handle); end
if ishandle(ad.mean_shape_handle); delete(ad.mean_shape_handle); end;
if ishandle(ad.aligned_edge_handle); delete(ad.aligned_edge_handle); end;
if ishandle(ad.text_edge_handle); delete(ad.text_edge_handle); end;
if ishandle(ad.aligned_text_edge_handle); delete(ad.aligned_text_edge_handle); end;
if ishandle(ad.mean_text_edge_handle); delete(ad.mean_text_edge_handle); end;
s = cell(size(ad.edge, 1),1);
for i=1:size(ad.edge,1)
   s{i} = num2str(i); 
end
ad.edge_handle           =  plot(ad.handles.imageAxis, ad.edge(:,2), ad.edge(:,1), '-o', 'LineWidth', 1.5, 'Color', 'r', 'MarkerFaceColor', 'g', 'MarkerEdgeColor', 'r', 'MarkerSize', 5, 'HitTest', 'off');
%ad.edge_handle(end+1)    = plot(ad.edge(ad.marker_indices,2), ad.edge(ad.marker_indices,1), 'd', 'LineWidth', 2, 'MarkerSize', 10, 'MarkerEdgeColor', 'y', 'MarkerFaceColor', 'g');
ad.text_edge_handle      = text(ad.edge(:,2), ad.edge(:,1), s, 'Parent', ad.handles.imageAxis, 'Background', 'w', 'HitTest', 'off');

ad.aligned_edge_handle        =  plot(ad.handles.imageAxis, ad.aligned_edge(:,2), ad.aligned_edge(:,1), '-o', 'LineWidth', 1.5, 'Color', 'r', 'MarkerFaceColor', 'g', 'MarkerEdgeColor', 'r', 'MarkerSize', 5, 'HitTest', 'off');
%ad.aligned_edge_handle(end+1) = plot(ad.aligned_edge(ad.marker_indices,2), ad.aligned_edge(ad.marker_indices,1), 'd',  'LineWidth', 2, 'MarkerSize', 10, 'MarkerEdgeColor', 'y', 'MarkerFaceColor', 'g');
ad.aligned_text_edge_handle   = text(ad.aligned_edge(:,2), ad.aligned_edge(:,1), s, 'Parent', ad.handles.imageAxis, 'Background', 'w', 'HitTest', 'off');

ad.mean_shape_handle        =  plot(ad.handles.imageAxis, ad.mean_shape(:,2), ad.mean_shape(:,1), '-o','LineWidth', 1.5, 'Color', 'g', 'MarkerFaceColor', 'r', 'MarkerEdgeColor', 'g', 'MarkerSize', 5, 'HitTest', 'off');
%ad.mean_shape_handle(end+1) = plot(ad.mean_shape(ad.marker_indices,2), ad.mean_shape(ad.marker_indices,1), 'd',  'LineWidth', 2, 'MarkerSize', 10, 'MarkerEdgeColor', 'r', 'MarkerFaceColor', 'w');
ad.mean_text_edge_handle    = text(ad.mean_shape(:,2), ad.mean_shape(:,1), s, 'Parent', ad.handles.imageAxis, 'Background', 'w', 'HitTest', 'off');

setappdata(0,'SectorWarpingToolData',ad);
doToggleEdge(ad.handles.ShowEdgeChk);
return;

%%%%
%
%
%%%%
function deleteHandles(ad)
if ishandle(ad.sectors_handle)
    delete(ad.sectors_handle);
end
if ishandle(ad.aligned_sectors_handle)
    delete(ad.aligned_sectors_handle);
end
if ishandle(ad.warped_sectors_handle)
    delete(ad.warped_sectors_handle);
end
if ishandle(ad.control_marker_handle)
    delete(ad.control_marker_handle);
end
return;

%%%%
%
%
%%%%
function displaySectors()
ad =getappdata(0, 'SectorWarpingToolData');
deleteHandles(ad);
if ~isfield(ad, 'sector')
    return;
end
hold(ad.handles.imageAxis, 'on');
ad.sectors_handle        =  plot(ad.handles.imageAxis, ad.sector.y, ad.sector.x, '.', 'HitTest', 'off');
ad.aligned_sectors_handle=  plot(ad.handles.imageAxis, ad.aligned_sector.y, ad.aligned_sector.x, '.', 'HitTest', 'off');
ad.warped_sectors_handle =  plot(ad.handles.imageAxis, ad.warped_sector.y, ad.warped_sector.x, '.', 'HitTest', 'off');
setappdata(0,'SectorWarpingToolData',ad);
doToggleSectors(ad.handles.ShowSectorsChk);
%%%
%
%%%
function clickedInWindow(win, evd)
ad =getappdata(0, 'SectorWarpingToolData');
p = get(ad.handles.imageAxis, 'CurrentPoint');
p = p(1,1:2);
if(get(ad.handles.warped_pos_radio, 'Value'))
    [ind, c, val] = findClosestSector(ad.warped_sector_info, p);
else
    [ind, c, val] = findClosestSector(ad.sector_info, p);
end
ph = plot(ad.handles.imageAxis, c(1), c(2), 'o', 'MarkerSize', 10, 'MarkerFaceColor', 'y', 'MarkerEdgeColor', 'r');

orig_area = ad.sector_info(ind).Area;
warped_area = ad.warped_sector_info(ind).Area;
%These function calls below return major and minor ellipse axis info
%note: these ellipse axis info are measured from centroid to edge of
%ellipse. Therefore we need to *2 to get the length and width of ellipse.
[major, minor, theta] = gtlib_growthTensor2Params(ad.sector_info(ind).growthTensor);
[warped_major, warped_minor, warped_theta] = gtlib_growthTensor2Params(ad.warped_sector_info(ind).growthTensor);
% major = major*2;
% minor = minor*2;
% warped_major = warped_major*2;
% warped_minor = warped_minor*2;

% fprintf('\n\n**********************\n');
% fprintf('Original Area: %6.6f\nWarped Area: %6.6f\nMajor Axis: %6.6f\nMinor Axis: %6.6f\nTheta: %6.6f\n', orig_area, warped_area, major, minor, theta);
% fprintf('**********************\n');
fprintf('\n\n**********************\n');
fprintf('Original Area: %6.6f\n',orig_area);
fprintf('Original Major Axis: %6.6f\n',major);
fprintf('Original Minor Axis: %6.6f\n',minor);
fprintf('Original Theta: %6.6f\n',theta);
fprintf('----------------------\n');
fprintf('Warped Area: %6.6f\n',warped_area);
fprintf('Warped Major Axis: %6.6f\n',warped_major);
fprintf('Warped Minor Axis: %6.6f\n',warped_minor);
fprintf('Warped Theta: %6.6f\n',warped_theta);
fprintf('**********************\n\n');
pause(1);
delete(ph);
setappdata(0, 'SectorWarpingToolData', ad);
return
%%%
%
%%%
function [ind, c, val] = findClosestSector(sector_info, p)
c = [sector_info.Centroid]';
c = reshape(c, 2, length(c)/2);
d = sqrt(sum((c - p(:)*ones(1, size(c,2))).^2,1));
[val, ind] = min(d);
val = val(1);
ind = ind(1);
c = c(:, ind);
return