function pp3dspot(immat,centre,width,pixscale,imcut,imcent)
% Function PP3DSPOT
%
%    pp3dspot(immat,centre,width,pixscale,imcut,imcent)
%
%    PP3DSPOT plots 3D images of individual pepperpot spots.  The function
%    requires a single pepperpot image as an input, extracts the data
%    surrounding the specified point, plots that data as a 3D coloured
%    contour plot and marks the centre position using a stem.
%
%    PP3DSPOT(IMMAT,CENTRE) extracts image data from the image file IMMAT,
%    surrounding the pixel location CENTRE, and plots that data.  IMMAT must
%    be a 2D matrix of pixel intensities; CENTRE a central pixel location
%    in [X Y].  The default width of the extracted image is 51 pixels.
%
%    PP3DSPOT(IMMAT,CENTRE,WIDTH) specifies the width of the sub-image that
%    is extracted.  WIDTH should be an odd number, since the sub-image is
%    specified in terms of a number of pixels either side of the central
%    pixel.  The default value is 51.
%
%    PP3DSPOT(IMMAT,CENTRE,WIDTH,PIXSCALE) gives a scale factor, PIXSCALE,
%    in [pixels-per-mm], that scales the axes correctly for the 3D plot.
%
%    PP3DSPOT(IMMAT,CENTRE,WIDTH,PIXSCALE,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.
%
%    PP3DSPOT(IMMAT,CENTRE,WIDTH,PIXSCALE,IMCUT,IMCENT) specifies the
%    actual pixel location, IMCENT, of the centre of the main image, in 
%    [X Y] pixels.
%

    if nargin < 2
        error('Need at least 2 input variables: pp3dspot(immat,centre,width,pixscale,imcut,imcent)') ;
    elseif nargin > 6
        error('Can only specify 6 input variables: pp3dspot(immat,centre,width,pixscale,imcut,imcent)') ;
    end

    if nargin < 6
        imcent = centre ;
    end
    if nargin < 5
        imcut = 0 ;
    end
    if nargin < 4
        pixscale = 1 ;
    end
    if nargin < 3
        width = 51 ;
    end

    spotrad = round( ( width - 1)./2 ) ;
    spotro = round( centre ) ;
    subim1 = immat(((spotro(2) - spotrad):(spotro(2) + spotrad)), ((spotro(1) - spotrad):(spotro(1) + spotrad))) ;

%    Apply cut to image intensity

    cutmat = ( immat > imcut ) ;
    immat = ( immat - imcut ).*cutmat ;
    imsize = size(immat) ;

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

    subim = immat(((spotro(2) - spotrad):(spotro(2) + spotrad)), ((spotro(1) - spotrad):(spotro(1) + spotrad))) ;
    spotrem = centre - spotro ;
    spotcent = ([spotrad + spotrem(1) + 1, spotrad + spotrem(2) + 1]) ;

    immax = max(max(subim1)) ;
    immin = min(min(subim1)) ;
    plotmax = ( round( (immax.*1.2)./100 ) ).*100 ;

%    Plot spot image

    minijet = colormap(jet(200)) ;
    minijet(1,:) = 1 ;

    figure, hold on ;
    surfc(flipud(subim)) ; %,'EdgeColor','interp') ;
    colormap(minijet) ;
    stem3(spotcent(1),spotcent(2),plotmax,'fill')
    axis tight ; axis square ;
    hold off ;

%    figure, imagesc(subim1,[immin immax]), axis image, hold on ;
%    colormap(gray) ;
%    plot(spotcent(1),spotcent(2),'xk') ;
%    hold off ;

%    Find image centre

%    cent = [(imsize(2) - 1)./2 (imsize(1) - 1)./2] ;





