function [im cr over]=ipe(Fx,Fy,x,y,sig,Lx,th,ybp,ytp,L,Nx,Ny)

Np=numel(x);

[cr over]=pgrid(Ny-y/L*Nx,x/L*Nx,Ny,Nx,[1 Nx 1 Ny],Np,sig/L*Nx,0);
over(over==0)=Np+1;


F=[Fx(:,1:Np)-Fx(:,1:Np)' Fx(:,Np+1:end)]+i*[Fy(:,1:Np)-Fy(:,1:Np)' Fy(:,Np+1:end)];
%dd=mdiff(-[x 0 Lx 0 0 0 0]+i*[y 0 0 0 0 ybp ytp])/L*Ny;
dd=[mdiff(-x+i*y) -x'.*exp(-i*th) (Lx-x').*exp(i*th) zeros(Np,2) i*(y'-ytp) i*(y'-ybp)];

F=abs(F).*dd./abs(dd);
F(isnan(F))=0;

[Fs ii]=(sort(F,2,'Descend'));

%for np=1:Np; Fs(np,:)=angle(F(np,ii(np,:)));end

im=0;
for np=1:6;
  Fn=[Fs(:,np); 0];
  im=im+abs(Fn(over)).*LineForce(-cr.*exp(i*(pi/2+angle(Fn(over)))),1/2,.1);
end
