float rand1(int x, int y) {
    x=(x+y*57);
    x=(x<<13)^x;
    return (1.0f-((x*(x*x*15731+789221)+1376312589)&0x7fffffff)/1073741824.0f);
}

float interpolate(float a, float b, float x) {
    x=(1-cos(x*PI))*0.5f;
    return a*(1-x)+b*x;
}


float smooth(int x, int y) {
    float corn = (rand1(x-1,y-1)+rand1(x+1,y-1)+rand1(x+1,y+1)+rand1(x-1,y+1))/16.0f;
    float side = (rand1(x-1,y)+rand1(x+1,y)+rand1(x,y-1)+rand1(x,y+1))/8.0f;
    return corn+side+rand1(x,y)/4.0f;
}

float sAndI(float x, float y) {
    int integer_X  = (int)x;
    float fractional_X = x - integer_X;
    int integer_Y = (int)y;
    float fractional_Y = y - integer_Y;
    float v1,v2,v3,v4,i1,i2;
    v1 = smooth(integer_X,     integer_Y);
    v2 = smooth(integer_X + 1, integer_Y);
    v3 = smooth(integer_X,     integer_Y + 1);
    v4 = smooth(integer_X + 1, integer_Y + 1);
    i1 = interpolate(v1 , v2 , fractional_X);
    i2 = interpolate(v3 , v4 , fractional_X);
    return interpolate(i1 , i2 , fractional_Y);
}

float noise(float persistance, int iter, int x, int y) {
    float tot,freq,amp;
    tot=0; freq=amp=1;
    for (int c=0;c<iter;c++) {
        tot+=sAndI(x*freq,y*freq)*amp;
        freq*=2; amp*=persistance;
    };
    return tot;
};