--- original_ppfit.m	2012-03-22 09:18:22.000000000 +0100
+++ ppfit.m	2012-03-22 16:17:53.000000000 +0100
@@ -113,17 +113,13 @@
            "Input data value pairs must have the same size.")
   endif
   ## Use linear least squares for the initial guess
-  a = zeros (numel (x), numel (xi));
-  for ii = 1:numel (xi)
-    yi = zeros (numel (xi), 1);
-    yi(ii) = 1;
-    a(:, ii) = interp1 (xi, yi, x, method, "extrap");
-  endfor
+  yi = eye (numel (xi));
+  a = interp1 (xi, yi, x(:), method, "extrap");
   ## y = (q*r) * yi + errors
   ws = warning ("off", "all"); 
   [q, r, p] = qr (a .* w, 0);
   [~, p] = sort (p);
-  yi = (r \ (q' * (y .* w))).'(p);
+  yi = (r \ (q' * (y .* w)))(p);
   warning (ws);
   if (any (strcmp (method, {"pchip", "cubic"})))
     ## Use fsolve() to minimize the weighted LSE