function [outdat,outnames] = ppcalimdat(immat,cent,lpos,rpos,lhm,rhm)
% Function PPCALIMDAT
%
%    [angs,angname] = ppcalimdat(immat,lpos,rpos,cent,lhm,rhm)
%
%    PPCALIMDAT calculates the angles between various points of Pepperpot
%    data images for a particular data set, and outputs this data as a
%    single matrix
%
%    [OUTDAT,OUTNAMES] = PPCALIMDAT(IMMAT) loads the image matrix IMMAT and
%    displays it on screen: the user can then find the left and right hand
%    edges of the horizontal ruler, the centre of the image and the left
%    and right hand ends of the midlineby hand and input this information
%    at the prompt in the form of a matrix:
%
%       > Please enter left-hand ruler coordinates:   [113 1884]
%
%      Various properties of the image are then calculated and output in
%    the OUTDAT cell array.  The corresponding name of the entry in each
%    cell of the OUTDAT array is contained in the OUTNAMES cell array; this
%    does not have to be specified.
%
%    The image data should be preloaded using the PEPPERPLOTFUN function.
%
%    [OUTDAT,OUTNAMES] = PPCALIMDAT(IMMAT,[CENT]) ...
%
%    [OUTDAT,OUTNAMES] = PPCALIMDAT(IMMAT,[CENT],[LPOS],[RPOS]) ... and
%
%    [OUTDAT,OUTNAMES] = PPCALIMDAT(IMMAT,[CENT],[LPOS],[RPOS],[LHM],[RHM])
%    makes it possible to specify the pixel locations of the image centre,
%    the left- and right-hand ends of the ruler and the left- and
%    right-hand points of the midline, where, for example:
%
%        LPOS = [LX LY] ;    RPOS = [RX RY] ;
%
%      The coordinate positions of BOTH ends of the ruler and BOTH ends of
%    the midline must be specified otherwise the function will not run.
%
%    The function requires the following input variables:
%
%    IMMAT - The data/image file with the image information, including
%              the file extension and the path (absolute or relative) -
%              must be a string.
%

    if nargin < 1
        error('Need at least 1 input variable: [outdat,outnames] = ppcalimdat(immat,cent,lpos,rpos,lhm,rhm)') ;
    elseif nargin > 6
        error('Can only specify 6 input variables: [outdat,outnames] = ppcalimdat(immat,cent,lpos,rpos,lhm,rhm)') ;
    elseif ( ( nargin == 3 ) | ( nargin == 5 ) )
        error('Must specify coordinates in pairs for ruler or midline: [outdat,outnames] = ppcalimdat(immat,cent,lpos,rpos,lhm,rhm)') ;
    elseif nargout < 1
        error('Must specify output matrix: [outdat,outnames] = ppcalimdat(immat,cent,lpos,rpos,lhm,rhm)') ;
    elseif nargout > 2
        error('Can only specify 2 output variables: [outdat,outnames] = ppcalimdat(immat,cent,lpos,rpos,lhm,rhm)') ;
    elseif nargin == 4
        if ( length(lpos) ~= 2 ) | ( length(rpos) ~= 2 )
            error('Coordinates must contain X and Y pixels: lpos = [lx ly]') ;
        end
    elseif nargin == 6
        if ( length(lpos) ~= 2 ) | ( length(rpos) ~= 2 ) ( length(lhm) ~= 2 ) | ( length(rhm) ~= 2 )
            error('Coordinates must contain X and Y pixels: lpos = [lx ly]') ;
        end
    end

%    Set up matrix with output variable names

    outnames = {['cent']; ['lpos']; ['rpos']; ['lhm']; ['rhm']; ['rulen']; ...
        ['ruang']; ['lang']; ['rang']; ['cang']; ['mltilt']; ['rumlang']; ['mllen']} ;

%    Load image data for display and analysis

    immin = min(min(immat)) ; immax = max(max(immat)) ;
    gsim = figure ; imagesc(immat, [immin immax]), axis image ;
    colormap(gray) ;

    ppbw = edge(immat,'canny') ;
    bwim = figure ; imagesc(ppbw), axis image ;
    colormap(gray) ;

%    Ask for positions of various coordinates if not provided

    if nargin < 2

        disp(' ') ;
        cent = input('Please enter image centre coordinates:  ') ;
        disp(' ') ;

    end

    if nargin < 4

        disp(' ') ;
        lpos = input('Please enter left-hand ruler coordinates:  ') ;
        disp(' ') ;
        rpos = input('Please enter right-hand ruler coordinates:  ') ;
        disp(' ') ;

    end

    if nargin < 6

        disp(' ') ;
        lhm = input('Please enter left-hand midline coordinates:  ') ;
        disp(' ') ;
        rhm = input('Please enter right-hand midline coordinates:  ') ;
        disp(' ') ;

    end

%    Calculate image angles

    cldiff = cent - lpos ;   % Difference between left ruler end and central point
    crdiff = rpos - cent ;   % Difference between right ruler end and central point
    ruldiff = rpos - lpos ;  % Difference between end points of ruler
    rulgrad = ruldiff(2) ./ ruldiff(1) ;     % Gradient or ruler
    ruloff = lpos(2) - ( lpos(1).*ruldiff ) ;       % Offset of ruler
    rulmid = lpos + ruldiff./2 ;                    % Midpoint of ruler
    mcdiff = rulmid - cent ;     % Difference between ruler midpoint and centre

    rulen = sqrt( sum(ruldiff.^2) ) ;   % Length of ruler
    cllen = sqrt( sum(cldiff.^2) ) ;    % Distance between central point and left ruler edge
    crlen = sqrt( sum(crdiff.^2) ) ;    % Distance between central point and right ruler edge
    mclen = sqrt( sum(mcdiff.^2) ) ;    % Distance between central point and ruler midpoint
    ruang = atan2( ruldiff(2), ruldiff(1) ) ;   % Angle of ruler
    clang = atan2( cldiff(2), cldiff(1) ) ;     % Angle of Centre-to-ruler-left line
    crang = atan2( crdiff(2), crdiff(1) ) ;     % Angle of Centre-to-ruler-right line
    mcang = atan2( mcdiff(2), mcdiff(1) ) ;     % Angle of midpoint-to-centre line

    lang = clang - ruang ;             % Angle between left-hand triangle edge and ruler
    rang = crang - ruang ;             % Angle between right-hand triangle edge and ruler
    cang = (pi - abs(lang) - abs(rang) ) ;  % Angle between left- and right-hand triangle edges
    mang = mcang + ruang ;                  % Angle between ruler and midpoint-to-centre line

    mllen = sqrt( (rhm(2) - lhm(2)).^2 + (rhm(1) - lhm(1)).^2 ) ;   % Length of midline
    mltilt = atan2( (rhm(2) - lhm(2)),(rhm(1) - lhm(1)) ) ;         % Tilt angle of midline
    rumlang = mltilt - ruang ;

%    Output data matrix

    outdat = {cent; lpos; rpos; lhm; rhm; rulen; ruang; lang; rang; ...
        cang; mltilt; rumlang; mllen} ;

    figure(gsim), hold on ;
    plot([lpos(1),rpos(1)], [lpos(2),rpos(2)],'-xr') ;
    plot([lpos(1),cent(1)], [lpos(2),cent(2)],'-xr') ;
    plot([cent(1),rpos(1)], [cent(2),rpos(2)],'-xr') ;
    plot([lhm(1),rhm(1)], [lhm(2),rhm(2)],'-xr') ;
    plot([rulmid(1),cent(1)], [rulmid(2),cent(2)],'-xc') ;
    hold off ;

    return








