
struct VsOutput
{
  float4 pos: POSITION;
  float4 tcp: TEXCOORD0;
  float4 tcn: TEXCOORD1;
  float2 ntc: TEXCOORD2;
};

float4 texsz : register(c16);
float4 pixofs: register(c17);
float4 ntmu  : register(c18);
float4 ntmv  : register(c19);

VsOutput vsh(float2 pos : POSITION)
{
  VsOutput O;
  float2 p=pos;
  O.pos=float4(p, 1, 1);
  float2 t=p*texsz.xy+texsz.zw;
  O.tcp=t.xyxy+float4(0, 0, pixofs.xy);
  O.tcn=t.xyxy-float4(0, 0, pixofs.xy);
  O.ntc.x=dot(float3(t, 1), ntmu.xyz);
  O.ntc.y=dot(float3(t, 1), ntmv.xyz);
  return O;
}


sampler    flowMap : register(s0);
sampler erosionMap : register(s1);
sampler   noiseMap : register(s2);

float4    pixelOfs : register(c0);
float4   params[2] : register(c1);


void psh(VsOutput I, out float4 newFlow : COLOR0, out float4 newErosion : COLOR1)
{
  float4 f=tex2D(   flowMap, I.tcp.xy);
  float4 e=tex2D(erosionMap, I.tcp.xy);
  float4 noise=tex2D(noiseMap, I.ntc.xy);
  float v=f.r, s=f.g, w=f.b, h=f.a;

  float Kq=params[0].x;
  float Kw=params[0].y;
  float Kr=params[0].z;
  float Kd=params[0].w;
  float minSlope=params[1].x;
  float       Kg=params[1].y;
  float       Kp=params[1].z;

  // calc gradient dir
  float hxp=tex2D(flowMap, I.tcp.zy).a;
  float hxn=tex2D(flowMap, I.tcn.zy).a;
  float hzp=tex2D(flowMap, I.tcp.xw).a;
  float hzn=tex2D(flowMap, I.tcn.xw).a;
  float dx=(hxp>hxn ? hxp-h : h-hxn);
  float dz=(hzp>hzn ? hzp-h : h-hzn);

  float2 gdir=normalize(float2(dx, dz));

  // get "bottom" height
  float bh=tex2D(flowMap, I.tcp.xy-gdir*pixelOfs.xy).a;

  #if !ERODE
  // propagate "top" data here
  float4 tf=tex2D(flowMap, I.tcp.xy+gdir*pixelOfs.xy);

  tf.a-=h;
  tf =(tf.a>=0 ? tf : float4(0, s, w, max(h-bh, minSlope)));
  float tdh=tf.a;

  float rain=(noise.x<=Kp)?noise.x:0;
  //float rain=Kp*0.5;

  v=sqrt(tf.r*tf.r+Kg*tdh);
  s=tf.g;
  w=tf.b*Kw+rain;
  #endif

  #if !FLOW
  // deposit/erode
  float dh=h-bh;
  float ds;

  if (dh<0)
  {
    ds=0.001-dh;
    if (ds>s) { ds=s; w=0; }
    v=0;
  }
  else
  {
    // calc transport capacity
    float slope=dh;
    //float slope=dh/sqrtf(dh*dh+1);
    float q=max(slope, minSlope)*v*w*Kq;

    // deposit/erode
    ds=s-q;
    ds=(ds>=0 ? ds*Kd : max(ds*Kr, -dh));
  }

  s-=ds;
  h+=ds;

  #if 1
    e.y+=ds;
    e.x-=min(e.y, 0);
    e.y =max(e.y, 0);
  #else
    e.y+=max(ds, 0);
    e.x-=min(ds, 0);
  #endif

  #endif // !FLOW

  // store results
  f.r=v; f.g=s; f.b=w; f.a=h;
  newFlow   =f;
  newErosion=e;
}
