// Ported from http://cs.brown.edu/people/pfelzens/dt/

sclass DistanceTransform {

  static float INF = 1E20f;

  // dt of 1d function using squared distance

  float[] dt(float[] f, int n) {

    float[] d = new float[n], z = new float[n+1];

    int[] v = new int[n];

    int k = 0;

    z[0] = -INF;

    z[1] = INF;

    for (int q = 1; q <= n-1; q++) {

      float s = ((f[q]+sqr(q))-(f[v[k]]+sqr(v[k])))/(2*q-2*v[k]);

      while (s <= z[k]) {

        k--;

        s = ((f[q]+sqr(q))-(f[v[k]]+sqr(v[k])))/(2*q-2*v[k]);

      }

      k++;

      v[k] = q;

      z[k] = s;

      z[k+1] = INF;

    }

    k = 0;

    for (int q = 0; q <= n-1; q++) {

      while (z[k+1] < q)

        k++;

      d[q] = sqr(q-v[k]) + f[v[k]];

    }

    ret d;

  }

  // dt of 2d function using squared distance

  void dt(FloatBWImage im) {

    int width = im.getWidth();

    int height = im.getHeight();

    float[] f = new[max(width,height)];

    // transform along columns

    for x to width: {

      for y to height:

        f[y] = im.getPixel(x, y);

      float[] d = dt(f, height);

      for y to height:

        im.setPixel(x, y, d[y]);

    }

    // transform along rows

    for y to height: {

      for x to width:

        f[x] = im.getPixel(x, y);

      float[] d = dt(f, width);

      for x to width:

        im.setPixel(x, y, d[x]);

    }

  }

  /* dt of binary image using squared distance */

  FloatBWImage dt(BWImage im, float threshold default .5f, bool brightIsOn) {

    int width = im.getWidth();

    int height = im.getHeight();

    var out = new FloatBWImage(width, height);

    for y to height:

      for x to width:

        out.setPixel(x, y,

          (im.getFloatPixel(x, y) >= threshold == brightIsOn)

            ? 0 : INF);

    dt(out);

    ret out;

  }

}