function [a] = N(m,mp,w)
pp = size(w);
for i=1:pp(2)
 if ((w(i)==m)|(w(i)==mp)) 
  if (m==mp)
   a(i) = pi^2;
  else
   a(i) = 0;
  end;
 else
   a(i) = (sin(pi*w(i))^2)/(w(i)-m)/(w(i)-mp);
 end;
end;
% normalize it to a delta-function
a = a / (pi^2);
end;

