scale = 2.0 / N;
for (j=-N/2;j<N/2;j++) {
   for (i=-N/2;i<N/2;i++) {
      p[0].y = j * scale;
      p[1].y = (j+1) * scale;
      p[2].y = p[1].y;
      p[3].y = p[0].y;
      p[0].x = i * scale;
      p[1].x = p[0].x;
      p[2].x = (i+1) * scale;
      p[3].x = p[2].x;
      for (k=0;k<4;k++)
         r[k] = p[k].x*p[k].x+p[k].y*p[k].y;
      if (r[0] > 1 || r[1] > 1 ||
          r[2] > 1 || r[3] > 1)
         continue;
      for (k=0;k<4;k++)
         p[k].z = sqrt(1 - r[k]);
      draw the polygon p[0], p[1], p[2], p[3]
   }
}