function varargout=wavelength_elastography(dpl,dx,dz,Pz,varargin)

    if nargin==4
        axe='Z';
    elseif nargin>4
        axe=varargin{1};
    end

    switch axe
        case 'Z'
            strain=zeros(size(dpl))*NaN;
            strain(1:end-Pz,:,:)=(dpl(1+Pz:end,:,:)-dpl(1:end-Pz,:,:))/(Pz*dz);
            autocorrStrain = sum(strain.^2,3)/size(strain,3);%5*

        case 'X'
            dpl=permute(dpl,[2 1 3]);
            strain=zeros(size(dpl))*NaN;
            strain(1:end-Pz,:,:)=(dpl(1+Pz:end,:,:)-dpl(1:end-Pz,:,:))/(Pz*dz);
            dpl=permute(dpl,[2 1 3]);
            strain=permute(strain,[2 1 3]);
            autocorrStrain=sum(strain.^2,3)/size(strain,3);%5/2*

        case 'XZ'
            strain=zeros(size(dpl))*NaN;
            strain(1:end-Pz,:,:)=(dpl(1+Pz:end,:,:)-dpl(1:end-Pz,:,:))/(Pz*dz);
            autocorrStrainZ = sum(strain.^2,3)/size(strain,3);
            dpl=permute(dpl,[2 1 3]);
            strain=zeros(size(dpl))*NaN;
            strain(1:end-Pz,:,:)=(dpl(1+Pz:end,:,:)-dpl(1:end-Pz,:,:))/(Pz*dz);
            dpl=permute(dpl,[2 1 3]);
            strain=permute(strain,[2 1 3]);
            autocorrStrainX=sum(strain.^2,3)/size(strain,3);
            autocorrStrain=autocorrStrainZ+autocorrStrainX; 
    end

    autocorrDPL=sum(dpl.^2,3)/size(dpl,3);
    LambdaZ = 2*pi*sqrt(autocorrDPL./autocorrStrain);

    if nargout==0
        x=size(dpl,2)*dx;
        z=size(dpl,1)*dz;
        figure,imagesc((x-mean(x)),z,real(LambdaZ));
        colorbar;
        title('The Shear Wavelength Map ')
        title(colorbar,'[mm]','FontSize',20)
        xlabel('Width [mm]','FontSize',20); ylabel('Depth [mm]','FontSize',20);
    elseif nargout==1
        varargout{1}=real(LambdaZ);
    end
    
end