
function [handlelHA,wimg,net_offset,handlelHAP,statistics]=PDMWalk_ShowModel(...
    aam, model, oversize_by, handlelHA, pictureH, pictureflag,...
    pts, handlelHAP,img2,w_image_pts,w_image_itself, w_image_indexes,speedflag)
if nargin<3
    oversize_by=1.4;
end
if nargin<4
    handlelHA=[];
end
if nargin<6
    pictureflag=0;
    pts=[];
end
if nargin<8
    handlelHAP=0;
end
if nargin<9
    img2=[];
end
if nargin<10
    w_image_pts=[];
end
if nargin<11
    w_image_itself=[];
end
if nargin<12
    w_image_indexes=[];
end
if nargin<13
    speedflag=[];
end
bWA=aam.bWA;
img = double(aam.Am)+aam.P * bWA * aam.multi;
img = reshape(img, [aam.siz(2), aam.siz(1), 3]);
% now we need to warp it to the following shape
rmodel=reshape(model',2,length(model)/2)';
%
% set up the display image using the original mean pointmodel
% zero the offset  (Data.plot_pdm)
mean_rmodel=reshape(aam.plot_pdm.Xm',2,length(aam.plot_pdm.Xm)/2)';
%mean_rmodel=reshape(aam.pdm.Xm',2,length(aam.pdm.Xm)/2)';
x_zoffset=min(mean_rmodel(:,1));
y_zoffset=min(mean_rmodel(:,2));
mean_zmodel(:,1)=mean_rmodel(:,1)-x_zoffset;
mean_zmodel(:,2)=mean_rmodel(:,2)-y_zoffset;
% centre into a larger display image
%oversize_by=1.2;
disp_img=zeros(floor([oversize_by*max(mean_zmodel(:,2)),oversize_by*max(mean_zmodel(:,1)),3]));
[rows,cols,planes]=size(disp_img);

x_zoffset=min(rmodel(:,1));
y_zoffset=min(rmodel(:,2));
zmodel(:,1)=rmodel(:,1)-x_zoffset;
zmodel(:,2)=rmodel(:,2)-y_zoffset;
x_middle_zmodel=mean(zmodel(:,1));
y_middle_zmodel=mean(zmodel(:,2));
x_offset=floor((cols/2)-(x_middle_zmodel)); %floor((cols-max(mean_zmodel(:,1)))/2);
y_offset=floor((rows/2)-(y_middle_zmodel)); %floor((rows-max(mean_zmodel(:,2)))/2);
net_offset=[x_offset-x_zoffset, y_offset-y_zoffset];

% Find the pixels that will be warped
timg=triwarp(100*ones(size(img(:,:,1))),aam.ompts,zmodel);
ind=find(timg>1);
% and offset the indexes into the display image
[indy,indx]=ind2sub(size(timg(:,:,1)),ind);
indx_off=indx+x_offset;
indx_off=min(indx_off,cols);
indx_off=max(indx_off,1);
indy_off=indy+y_offset;
indy_off=min(indy_off,rows);
indy_off=max(indy_off,1);
% convert from subscripts back to array index
ind_off=sub2ind(size(disp_img(:,:,1)),indy_off,indx_off);
% now warp the mean image (appearance) to right size
wimg=triwarp(img,aam.ompts,zmodel); % upsample to zero offset
% and insert the pixels that have been warped into the display
[wrows,wcols,wplanes]=size(wimg);
rdisp_img=reshape(disp_img,rows*cols,planes);
rwimg=reshape(wimg,wrows*wcols,wplanes);
rdisp_img(ind_off,:)=rwimg(ind,:);
disp_img=reshape(rdisp_img,[rows,cols,planes]);
% insert into the display image
wimg = disp_img;
wimg = double(disp_img)/255;
%wimg = wimg - min(min(min(wimg)));  % normalise the image brightness
%wimg = wimg / max(max(max(wimg)));
if ~isempty(img2)
    if speedflag==1 % only compute statistics using middle of warped original image
        statistics.n_vertices=length(model(:));
        statistics.dist_pts=sqrt(sum((model(:)-w_image_pts(:)).^2))/statistics.n_vertices;
        [rows,cols,planes]=size(w_image_itself);
        rimg3=reshape(w_image_itself,rows*cols,planes);
        rimg=reshape(img,rows*cols,planes);
        statistics.n_pixels=planes*length(w_image_indexes);
        statistics.dist_pts=sqrt(sum((rimg3(w_image_indexes,:)-rimg(w_image_indexes,:)).^2))/statistics.n_pixels;
    else % put the newly generated image into the original picutre
        %M2 = calctrm(pts, model, ones(length(pts)/2, 1));
        % the model comes from PDMWalk_PlotModel and accounts for slider positions
        % It will use scaling,rotation,translation (even though not all
        % might have been incorporated into the model). Otherwise one would
        % need to use pmalign.m
        M2 = calctrm(pts, model(:), ones(length(pts)/2, 1));
        model2=applytrm(model(:), M2);
        statistics.n_vertices=length(model2(:));
        statistics.dist_pts=sqrt(sum((model2(:)-pts(:)).^2))/statistics.n_vertices;
        picture_rmodel=reshape(model2',2,length(model)/2)';
        % ensure the image fills the full range of intensities
        img2=double(img2); % imported image
        %img2=img2-min(min(min(img2))); %*****
        %img2=255*img2/max(max(max(img2)));%*****
        img=double(img); % mean or reference image
        %img=img-min(min(min(img)));%*****
        %img=255*img/max(max(max(img)));%*****
        % Find the pixels that will be warped
        timg=triwarp(100*ones(size(img2(:,:,1))),aam.ompts,picture_rmodel);
        [rows2,cols2,planes2]=size(img2);
        timg=imerode(timg,strel('disk',7)); %shrink the area a little to prevent black margins *****
        ind=find(timg>1);
        [ii,jj]=ind2sub(size(timg),ind); %correct indexes but ind is wrong size
        % now warp the mean image (appearance) to right size
        picture_img=triwarp(img,aam.ompts,picture_rmodel); % upsample to zero offset
        % and insert the pixels that have been warped into the display
        [prows,pcols,pplanes]=size(picture_img);
        temp=zeros(size(img2)); % placing the model in the same size image
        temp(1:min(rows2,prows),1:min(cols2,pcols),:)=picture_img(1:min(rows2,prows),1:min(cols2,pcols),:);
        ind=find(ii>rows2);
        ii(ind)=rows2;
        ind=find(jj>cols2);
        jj(ind)=cols2;
        % find the indexes for the correct indexes
        indtemp=sub2ind(size(temp),ii,jj);
        % now reshape so as to use the indexes in ind
        rtemp=reshape(temp,rows2*cols2,planes2);
        rimg2=reshape(img2,rows2*cols2,planes2);
        rimg2(indtemp,:)=rtemp(indtemp,:);
        picture_wimg=reshape(rimg2,[rows2,cols2,planes2]);
    end
else
    statistics=[];
    dist_pts=0;
    dist_aam=0;
    picture_wimg=[];
end
%if nargout==1

if isempty(handlelHA)
    handlelHA=imagesc(wimg);
else
    if ishandle(handlelHA)
        set(handlelHA,'CData',wimg);
    else
        handlelHA=image(wimg);
    end
end

axis equal
%drawnow
if nargin>3
    if get(pictureH(2),'value') | get(pictureH(5),'value')
        set(handlelHA,'visible','on')
    else
        set(handlelHA,'visible','off')
    end
end
if pictureflag~=0
    hold on
    if isempty(handlelHAP)
        handlelHAP=image(uint8(picture_wimg));
    else
        if ishandle(handlelHAP)
            set(handlelHAP,'CData',uint8(picture_wimg));
        else
            handlelHAP=image(uint8(picture_wimg));
        end
    end
    hold off
    axis equal
    %drawnow
    if nargin>3
        value=get(pictureH(3),'value');
        if value<1
            set(handlelHAP,'visible','off')
        else
            set(handlelHAP,'visible','on')
        end
    end
end
