function [y1,y2,x] = strehl1(p,w)

s = size(p);
[z zmax] = psf(p,zeros(s));
ncen = s(1)/2+1;
M = 100;
x = zeros(M+1,1);
y1 = x;
y2 = x;

for k=1:M+1
    x(k)=3*(k-1)/M;
    [z zs] = psf(p,x(k)*w,zmax);
    y1(k)= zs;
    y2(k)=z(ncen,ncen);
end