import java.applet.*; import java.awt.*; import java.awt.image.*; import java.awt.event.*; import java.io.*; import java.net.*; import java.text.*; import java.util.*; import java.util.zip.*; import netscape.javascript.*; import javax.comm.*; import javax.sound.midi.*; import javax.sound.midi.spi.*; import javax.sound.sampled.*; import javax.sound.sampled.spi.*; import javax.xml.parsers.*; import javax.xml.transform.*; import javax.xml.transform.dom.*; import javax.xml.transform.sax.*; import javax.xml.transform.stream.*; import org.xml.sax.*; import org.xml.sax.ext.*; import org.xml.sax.helpers.*; public class MapImage1 extends BApplet {PSystem ps;
Arcball arcball;

ImageMap iMap;
float[]  input = new float[2];
float[]  output = new float[1];

SteepestMap sMap;

void setup() {

  size(400, 600);

  iMap = new ImageMap(1, "small.gif");
  iMap.interpMode = INTERP_BICUBIC;
  
  sMap = new SteepestMap("small.gif");
  sMap.interpMode = INTERP_BICUBIC;
}

void loop() {

  background(255);

  for (int x = 0; x < width; x++) {
    for (int y = 0; y < height; y++) {

      input[0] = (float)x / width;
      input[1] = (float)y / height;

      iMap.map(input, output);

      stroke((int)(output[0] * 255));
      point(x, y);

    }

  }

  stroke(255, 0, 0, 100);
  sMap.draw2DVectors(width, height, 100, 100, 20);

}

class SinCosField2D extends Map {

  SinCosField2D() {
    super(2, 2);
  }

  void map(float[] in, float[] out) {
    out[0] = 0.5f - ((sin(TWO_PI*(in[0] + (float)mouseX / width)) + 1) / 2) * ((cos(frame / 100.0f + TWO_PI*in[1]) + 1) / 2);
    out[1] = 0.5f - ((cos(frame / 77.0f + TWO_PI*in[0]) + 1) / 2) * ((sin(TWO_PI*(in[1] + (float)mouseY / height)) + 1) / 2);
    //out[0] = (0.5 - (float)mouseX / width) * in[0];
    //out[1] = (0.5 - (float)mouseY / height) * in[1];
  }
  
}

class SqField2D extends Map {

  SqField2D() {
    super(2, 2);
  }

  void map(float[] in, float[] out) {
    out[0] = sqrt(in[0]);
    out[1] = in[1] * in[1];
  }
  
}

class SinCosField3D extends Map {

  SinCosField3D() {
    super(3, 3);
  }

  void map(float[] in, float[] out) {
    out[0] = ((sin(TWO_PI*(in[0] + (float)mouseX / width)) + 1) / 2) * ((cos(frame / 100.0f + TWO_PI*in[1]) + 1) / 2);
    out[1] = ((cos(frame / 77.0f + TWO_PI*in[0]) + 1) / 2) * ((sin(TWO_PI*(in[1] + (float)mouseY / height)) + 1) / 2);
    out[2] = in[2] * in[2];
  }

}


public class MapForce extends Force {
    private PSystem pSystem;
    
    public Map mapObj;
    public float inputDivisor;
    public float outputMultiplier = 1.0f;

    MapForce() {
      super();
      inputDivisor = (width + height)/2.0f;
    }

    MapForce(Map map) {
      this();
      this.mapObj = map;
    }

  void init(PSystem pSystem) {
    this.pSystem = pSystem;
  }


    void applyForce() {
      Particle p;
      Particle particles[] = pSystem.particles;
      for (int i = 0; i < pSystem.nParticles; i++) {
        p = particles[i];
        applyForce(p);
      }
    }

    void applyForce(Particle p) {
       float[] mapOut = new float[3];
       float[] mapIn = new float[3];
       mapIn[0] = p.pos[0] / inputDivisor;
       mapIn[1] = p.pos[1] / inputDivisor;
       mapIn[2] = p.pos[2] / inputDivisor;
       
       mapObj.map(mapIn, mapOut);
       
       mapOut[0] = mapOut[0] * outputMultiplier;
       mapOut[1] = mapOut[1] * outputMultiplier;
       mapOut[2] = mapOut[2] * outputMultiplier;
       
       p.addForce(mapOut[0], mapOut[1], mapOut[2]);
    }
}


final static int INTERP_NONE = 0;
final static int INTERP_BILINEAR = 1;
final static int INTERP_BICUBIC = 2;

class ImageMap extends Map {

  BImage image;
  int imageWidth, imageHeight;
  float xStep, yStep;
  int   interpMode = INTERP_BICUBIC;

  // B-Spline blending
  final float[][] ms_aafBSBlend = {
    { 1.0f/6.0f, -3.0f/6.0f,  3.0f/6.0f, -1.0f/6.0f },
    { 4.0f/6.0f,  0.0f/6.0f, -6.0f/6.0f,  3.0f/6.0f },
    { 1.0f/6.0f,  3.0f/6.0f,  3.0f/6.0f, -3.0f/6.0f },
    { 0.0f/6.0f,  0.0f/6.0f,  0.0f/6.0f,  1.0f/6.0f }
  };

  // Catmull-Rom blending (doesn't quite work)
  final float[][] ms_aafCRBlend = {
    { 0.0f, -0.5f,  1.0f, -0.5f },
    { 1.0f,  0.0f, -2.5f,  1.5f },
    { 0.0f,  0.5f,  2.0f, -1.5f },
    { 0.0f,  0.0f, -0.5f,  0.5f }
  };
  
  final float[][] m_aafBlend = ms_aafBSBlend;

  ImageMap(int outDim) {
    super(2, outDim);
  }

  ImageMap(int outDim, String imageName) {
    super(2, outDim);
    load(imageName);
  }

  void load(String imageName) {
    image = loadImage(imageName);
    imageWidth = image.width;
    imageHeight = image.height;
    xStep = 1.0f / imageWidth;
    yStep = 1.0f / imageHeight;
  }

  int xPixel(float x) {
    return max(0, min(imageWidth - 1, (int)(x * imageWidth)));
  }

  int yPixel(float y) {
    return max(0, min(imageHeight -1, (int)(y * imageHeight)));
  }

  float pixelValue(int x, int y) {
    return brightness(image.pixels[max(0, min(imageHeight-1,y)) * imageWidth + max(0, min(imageWidth-1,x))]);
  }
  
  float pixelValue(float x, float y) {
    return brightness(image.pixels[yPixel(y) * imageWidth + xPixel(x)]);
  }

  float bilinear(float x, float y) {
    int L = (int)(x / xStep);
    int K = (int)(y / yStep);

    float a = (x/xStep) - L;
    float b = (y/yStep) - K;

    int Lp1 = L + 1;
    if (Lp1 > imageWidth - 1) Lp1 = imageWidth - 1;

    int Kp1 = K + 1;
    if (Kp1 > imageHeight - 1) Kp1 = imageHeight - 1;

    return  (1-a)*(1-b)*pixelValue(L, K) + a*(1-b)*pixelValue(Lp1, K) + (1-a)*b*pixelValue(L, Kp1) + a*b*pixelValue(Lp1,Kp1);
  }

  // Bicubic algorithm taken from David Eberly's FANTASTIC www.magic-software.com site
  float bicubic(float fX, float fY) {
    // compute x-index and clamp to image
    float fXIndex = fX/xStep;
    int iX = (int)fXIndex;
    if ( iX < 0 || iX > imageWidth - 1 )
    return 0;

    // compute y-index and clamp to image
    float fYIndex = fY/yStep;
    int iY = (int)fYIndex;
    if ( iY < 0 || iY > imageHeight - 1 )
    return 0;

    float[] afU = new float[4];
    afU[0] = 1.0f;
    afU[1] = fXIndex - iX;
    afU[2] = afU[1]*afU[1];
    afU[3] = afU[1]*afU[2];

    float[] afV = new float[4];
    afV[0] = 1.0f;
    afV[1] = fYIndex - iY;
    afV[2] = afV[1]*afV[1];
    afV[3] = afV[1]*afV[2];

    // compute P = M*U and Q = M*V
    float[] afP = new float[4];
    float[] afQ = new float[4];
    int iRow, iCol;
    for (iRow = 0; iRow < 4; iRow++)
    {
      afP[iRow] = 0.0f;
      afQ[iRow] = 0.0f;
      for (iCol = 0; iCol < 4; iCol++)
      {
        afP[iRow] += m_aafBlend[iRow][iCol]*afU[iCol];
        afQ[iRow] += m_aafBlend[iRow][iCol]*afV[iCol];
      }
    }

    // compute (M*U)^t D (M*V) where D is the 4x4 subimage containing (x,y)
    iX--;
    iY--;
    float fResult = 0.0f;
    for (iRow = 0; iRow < 4; iRow++)
    {
      int iYClamp = iY + iRow;
      if ( iYClamp < 0 )
      iYClamp = 0;
      else if ( iYClamp > imageHeight - 1 )
      iYClamp = imageHeight - 1;

      for (iCol = 0; iCol < 4; iCol++)
      {
        int iXClamp = iX + iCol;
        if ( iXClamp < 0 )
        iXClamp = 0;
        else if ( iXClamp > imageWidth - 1 )
        iXClamp = imageWidth - 1;

        fResult += afQ[iRow]*pixelValue(iXClamp,iYClamp)*afP[iCol];
      }
    }

    return fResult;
  }
  
  float noInterp(float fX, float fY) {
    return pixelValue(fX, fY);
  }

  float valWithInterp(float fX, float fY) {
    if (interpMode == INTERP_BICUBIC) {
      return bicubic(fX, fY) / 255.0f;
    }
    else if (interpMode == INTERP_BILINEAR) {
      return bilinear(fX, fY) / 255.0f;
    }
    else {
      return noInterp(fX, fY) / 255.0f;
    }
  }

  // The default mapping just returns the brightness value
  void map(float[] in, float[] out) {
    out[0] = valWithInterp(in[0], in[1]);
  }
}

class SteepestMap extends ImageMap {

  int numDirTests = 32;
  float testStep;
  float stepDivisor = 8;

  SteepestMap() {
    super(2);
    isPolar = true;
  }

  SteepestMap(String imageName) {
    super(2);
    isPolar = true;
    load(imageName);
  }

  void load(String imageName) {
    super.load(imageName);
    testStep = (xStep + yStep) / 2;
  }

  void map(float[] in, float[] out) {
    float dirStep = TWO_PI / numDirTests;
    float bestDir = 0;
    float bestVal = 1.1f;
    float bestSize = 0;
    float fX = in[0];
    float fY = in[1];
    float currentBright = valWithInterp(fX, fY);
    
    float epsilon = testStep / stepDivisor;
    
    float testBright;
    
    for (float dir = 0; dir < TWO_PI; dir += dirStep) {
    
      testBright = valWithInterp(fX + epsilon * cos(dir), fY + epsilon * sin(dir));
      if (testBright < bestVal) {
        bestVal = testBright;
        bestDir = dir;
        bestSize = min(1, (currentBright - testBright) * stepDivisor);
      }    
    }
    
    if (bestSize < 0) {
      bestSize = 0;
    }
    
    out[0] = bestDir / TWO_PI;
    out[1] = bestSize;
  }
}

class Map {
  int inDim;
  int outDim;
  boolean isPolar = false;

  Map(int inDim, int outDim) {
    this.inDim = inDim;
    this.outDim = outDim;
  }

  void map(float[] in, float[] out) {
    for (int i = 0; i < outDim; i++) {
      out[i] = 0;
    }
  }

  void draw2DVectors(float fieldWidth, float fieldHeight, int xCells, int yCells, float baseLength) {
    if (isPolar) {
      draw2DVectorsPolar(fieldWidth, fieldHeight, xCells, yCells, baseLength);
    }
    else {
      draw2DVectorsCartesian(fieldWidth, fieldHeight, xCells, yCells, baseLength);
    }
  }

  void draw3DVectors(float fieldWidth, float fieldHeight, float fieldDepth, int xCells, int yCells, int zCells, float baseLength) {
    if (outDim != 3 || inDim != 3) {
      return;
    }

    float[] output = new float[3];
    float[] input = new float[3];

    float xInc = fieldWidth / xCells;
    float yInc = fieldHeight / yCells;
    float zInc = fieldDepth / zCells;

    float xIn, yIn, zIn;

    for (int x = 0; x < xCells; x++) {

      for (int y = 0; y < yCells; y++) {

        for (int z = 0; z < zCells; z++) {

          xIn = (x + 0.5f) / (float)xCells;
          yIn = (y + 0.5f) / (float)yCells;
          zIn = (z + 0.5f) / (float)zCells;
          
          input[0] = xIn;
          input[1] = yIn;
          input[2] = zIn;

          map(input, output);

          push();
          translate((x + 0.5f) * xInc, (y + 0.5f) * yInc, (z + 0.5f) * zInc);
          line(0, 0, 0, baseLength * output[0], baseLength * output[1], baseLength * output[2]);
          pop();
        }
      }
    }
  }

  void draw2DVectorsCartesian(float fieldWidth, float fieldHeight, int xCells, int yCells, float baseLength) {
    if (outDim != 2 || inDim != 2) {
      return;
    }

    float[] output = new float[2];
    float[] input = new float[2];

    float xInc = fieldWidth / xCells;
    float yInc = fieldHeight / yCells;

    float xIn, yIn;

    for (int x = 0; x < xCells; x++) {

      for (int y = 0; y < yCells; y++) {

        xIn = (x + 0.5f) / (float)xCells;
        yIn = (y + 0.5f) / (float)yCells;
        input[0] = xIn;
        input[1] = yIn;

        map(input, output);

        push();
        translate((x + 0.5f) * xInc, (y + 0.5f) * yInc);
        line(0, 0, baseLength * output[0], baseLength * output[1]);
        pop();
      }
    }
  }

  void draw2DVectorsPolar(float fieldWidth, float fieldHeight, int xCells, int yCells, float baseLength) {
    if (outDim != 2 || inDim != 2) {
      return;
    }

    float[] output = new float[2];
    float[] input = new float[2];

    float xInc = fieldWidth / xCells;
    float yInc = fieldHeight / yCells;

    float xIn, yIn;

    for (int x = 0; x < xCells; x++) {

      for (int y = 0; y < yCells; y++) {

        xIn = (x + 0.5f) / (float)xCells;
        yIn = (y + 0.5f) / (float)yCells;
        input[0] = xIn;
        input[1] = yIn;

        map(input, output);

        push();
        translate((x + 0.5f) * xInc, (y + 0.5f) * yInc);
        line(0, 0, output[1] * baseLength * cos(TWO_PI * output[0]), output[1] * baseLength * sin(TWO_PI * output[0]));
        pop();
      }
    }
  }

}


}