function [spotxydat] = ppspotallxy(immat,holepos,spotpos,width,imcut)
% Function PPSPOTALLXY.M
%
%    [spotxydat] = ppspotallxy(immat,holepos,spotpos,width,imcut)
%
%    PPSPOTALLXY calculates the X and Y offset, in [number of pixels], for
%    every single pixel within a given area from a specified reference
%    position.  This reference position should be the expected hole centre
%    in [X,Y].  The width of the search area, WIDTH, and an image intensity
%    low cut, IMCUT, must also be specified.
%
%    [SPOTXYDAT] = PPSPOTALLXY(IMMAT,HOLEPOS,SPOTPOS,WIDTH) finds the X and
%    Y pixel offsets, and the pixel intensities, for a given grayscale
%    image IMMAT.  The position of the origin from which the offset is to
%    be calculated is given by HOLEPOS, in [number of pixels].  The
%    position of the centre of the spot is given by SPOTPOS, also in
%    [number of pixels].  The size of the area in which to calculate is
%    given by WIDTH, also in [number of pixels].
%
%    The dimensions/arrangement of the matrices are as follows:
%
%        IMMAT - An [A x B] matrix of grayscale image intensity (usually
%                2048 x 2048).
%
%        HOLEPOS - The predicted hole position, [X Y], in [number of pixels].
%
%        SPOTPOS - The centre of the spot, [X Y], in [number of pixels].
%
%        WIDTH - The width of the area (in both X and Y) in which to
%                calculate the X, Y and intensity values, in [number of
%                pixels].  The default value is 51.
%
%        SPOTXYDAT - A [WIDTH x WIDTH x 3] matrix, where each entry (M,N,:)
%                    corresponds to the X, Y and Intensity data for a
%                    single pixel: [M,N,1:2] gives the horizontal and
%                    vertical pixel offset for each pixel from the
%                    specified hole position; [M,N,3] is the light
%                    intensity for that pixel.
%
%    [SPOTXYDAT] = PPSPOTALLXY(IMMAT,HOLEPOS,SPOTPOS,WIDTH,IMCUT) also
%    applies a cut to the image data, removing pixel intensities below the
%    intensity level IMCUT by subtracting the grayscale value IMCUT from
%    each pixel value, then making the minimum value of any pixel
%    zero.  The default value is IMCUT = 0.
%

    if nargin < 3
        error('Need at least 3 input variables: [spotxydat] = ppspotallxy(immat,holepos,spotpos,width,imcut)') ;
    elseif nargin > 5
        error('Can only specify 5 input variables: [spotxydat] = ppspotallxy(immat,holepos,spotpos,width,imcut)') ;
    elseif nargout ~= 1
        error('Must specify a single output variable: [spotxydat] = ppspotallxy(immat,holepos,spotpos,width,imcut)') ;
    end

    if nargin < 5
        imcut = 0 ;
    end
    if nargin < 4
        width = 41 ;
    end

%    holepos = round( holepos ) ;
    spotpos = round( spotpos ) ;
    width = round( width ) ;

    spotrad = round( ( width - 1)./2 ) ;

%    Apply cut to image intensity

    if imcut > 0
        cutmat = ( immat > imcut ) ;
        immat = ( immat - imcut ).*cutmat ;
    end

    imsize = size(immat) ;

    spotxydat = zeros(width,width,3) ;

    for i = 1:width
        witmat(:,i) = [1:1:width] ;
    end

    spotxydat(:,:,1) = witmat' ;
    spotxydat(:,:,2) = witmat ;

    offsets = spotpos - holepos ;
    spotxydat(:,:,1) = spotxydat(:,:,1) - (spotrad + 1) + offsets(1) ;
    spotxydat(:,:,2) = spotxydat(:,:,2) - (spotrad + 1) + offsets(2) ;

%    Select sub-pixel image from main matrix using spot position

    if ( ( spotpos(2) - spotrad ) <= 0 ) | ( ( spotpos(1) - spotrad ) <= 0 ) | ...
            ( ( spotpos(2) + spotrad ) > imsize(1) ) | ( ( spotpos(1) + spotrad ) > imsize(2) )

%        warning('Combination of SPOTPOS and WIDTH mean sub-image lies outside the image dimensions') ;

%        spotxydat(:,:,3) = 0 ;

    else

        subim = immat(((spotpos(2) - spotrad):(spotpos(2) + spotrad)), ((spotpos(1) - spotrad):(spotpos(1) + spotrad))) ;

        spotxydat(:,:,3) = subim ;

    end

    return