--- a/scripts/image/jet.m	Tue Apr 09 13:47:43 2013 -0400
+++ b/scripts/image/jet.m	Tue Apr 09 11:37:33 2013 -0700
@@ -45,17 +45,45 @@
 
   if (n == 1)
     map = [0, 1, 1];
-  elseif (n > 1)
-    x = linspace (0, 1, n)';
-    r = (x >= 3/8 & x < 5/8) .* (4 * x - 3/2) ...
-      + (x >= 5/8 & x < 7/8) ...
-      + (x >= 7/8) .* (-4 * x + 9/2);
-    g = (x >= 1/8 & x < 3/8) .* (4 * x - 1/2) ...
-      + (x >= 3/8 & x < 5/8) ...
-      + (x >= 5/8 & x < 7/8) .* (-4 * x + 7/2);
-    b = (x < 1/8) .* (4 * x + 1/2) ...
-      + (x >= 1/8 & x < 3/8) ...
-      + (x >= 3/8 & x < 5/8) .* (-4 * x + 5/2);
+  elseif (n == 2)
+    map = [0, 0, 1
+           0, 1, 1];
+  elseif (n > 2)
+    nel = ceil (n/4);           # number of elements
+    idx1 = ceil (3/8 * n) + 1;  # ~3/8*n for large n
+    if (mod (n, 8) == 2)
+      idx1++;
+    endif
+    idx2 = idx1 + nel - 1;      # ~5/8*n for large n
+    idx3 = min (idx2 + nel, n); # ~7/8*n for large n
+
+    r = zeros (n, 1);
+    r(idx1:idx2, 1) = [1:nel] / nel;
+    r(idx2:idx3, 1) = 1;
+    nel2 = n - idx3; 
+    r(idx3:(idx3+nel2), 1) = [nel:-1:(nel - nel2)] / nel;
+
+    idx1 = idx1 - nel;          # ~1/8*n for large n
+    idx2 = idx1 + nel - 1;      # ~3/8*n for large n
+    idx3 = min (idx2 + nel, n); # ~5/8*n for large n
+
+    g = zeros (n, 1);
+    g(idx1:idx2, 1) = [1:nel] / nel;
+    g(idx2:idx3, 1) = 1;
+    nel2 = min (nel, n - idx3); 
+    g(idx3:(idx3+nel2), 1) = [nel:-1:(nel - nel2)] / nel;
+
+    idx1 = max (idx2 - nel, 1); # ~1/8*n for large n
+    idx2 = idx2;                # ~3/8*n for large n
+    idx3 = idx3;                # ~5/8*n for large n
+
+    b = zeros (n, 1);
+    nel2 = min (nel, idx1-1); 
+    b(1:idx1, 1) = [(nel - nel2):nel] / nel;
+    b(idx1:idx2, 1) = 1;
+    nel2 = min (nel, n - idx3); 
+    b(idx2:(idx2+nel2), 1) = [nel:-1:(nel - nel2)] / nel;
+
     map = [r, g, b];
   else
     map = zeros (0, 3);