# Diamond Square Algorithm Implementation Problem

Im working on an implementation of the Diamon Square algortihm, and I can't quite seem to get it to look right.

I'm working off this tutorial http://www.lighthouse3d.com/opengl/terrain/index.php?mpd2 which might help explain some of my code.

The function that I know I need to change somehow is below.

It simply takes the min and max X/Z coordinates for the current square being subdivided (which index into an array of Y-Values where these min's and max's already have a Y value computed) and a displacement.

I'm pretty sure the step where I calculate a new midpoint is correct (setting E) but the part that I know is wrong is where I set the midpoint of each side (F G H I)

The forumla seems like it should work, but on high levels of recursion I see almost veritcal drops instead of a smoother transition:

http://img185.imageshack.us/img185/1323/failyh0.png

An alternate approach I tried that -almost- worked was just doing:

yMap[xMin][zMidpnt] = (A + C ) / 2;  // F
yMap[xMidpnt][zMin] = (A + B) / 2 ;  // G
yMap[xMax][zMidpnt] = (B + D)/ 2;  // H
yMap[xMidpnt][zMax] = (C + D) / 2 ;  // I

Which produces a more better landscape but isn't really performing the square step of the algorithm correctly (I think) because it does not add any random displacement to these midpoints.

Any help on how to fix this would be appreciated.

``````void genLandscapeRec(int xMin, int xMax, int zMin, int zMax, float disp) {

float nDisp =  (float)(disp * (float)pow(2, -(double)roughnes)); // New displacement

int xMidpnt = (xMin+xMax) / 2;
int zMidpnt = (zMin+zMax) / 2;

GLfloat A = yMap[xMin][zMin];
GLfloat B = yMap[xMax][zMin];
GLfloat C = yMap[xMin][zMax];
GLfloat D = yMap[xMax][zMax];

GLfloat E = (A + B + C + D) / 4 + randDisp(disp); // Midpoint height

yMap[xMidpnt][zMidpnt] = E; // Set E

yMap[xMin][zMidpnt] = (A + C + E) / 3 + randDisp(disp);  // F
yMap[xMidpnt][zMin] = (A + B + E) / 3 + randDisp(disp);  // G
yMap[xMax][zMidpnt] = (B + D + E) / 3 + randDisp(disp);  // H
yMap[xMidpnt][zMax] = (C + D + E) / 3 + randDisp(disp);  // I

if ((xMidpnt - xMin) > 1) { // Subdivide if new midpoint length will be > 1

genLandscapeRec(xMin, xMidpnt, zMin, zMidpnt, nDisp);
genLandscapeRec(xMidpnt, xMax, zMin, zMidpnt, nDisp);
genLandscapeRec(xMin, xMidpnt, zMidpnt, zMax, nDisp);
genLandscapeRec(xMidpnt, xMax, zMidpnt, zMax, nDisp);
}
}

float randDisp(float disp) {

return (((float)rand()) / ((float)RAND_MAX)) * (disp + (disp/3)) - disp/3;
}
``````
LVL 1
###### Who is Participating?
I wear a lot of hats...

"The solutions and answers provided on Experts Exchange have been extremely helpful to me over the last few years. I wear a lot of hats - Developer, Database Administrator, Help Desk, etc., so I know a lot of things but not a lot about one thing. Experts Exchange gives me answers from people who do know a lot about one thing, in a easy to use platform." -Todd S.

Author Commented:
The following globals and initialization function may also be helpful:
``````// fractal landscape globals
static GLfloat recDepth = 9;
static GLfloat roughnes = 1;
static GLfloat yMap[1025][1025];

void genLandscape(void) {

int dMax = (int)pow(2, recDepth);
int dMin = 0;

float disp = 3;

yMap[dMin][dMax] = 0.0;
yMap[dMin][dMax] = 0.0;
yMap[dMax][dMin] = 0.0;
yMap[dMax][dMax] = 0.0;

genLandscapeRec(dMin, dMax, dMin, dMax, disp);
}
``````
Commented:
hey iMr_Ki,

probably the randDisp-function is returning wrong (too high) values. check the values it is returning with an assert-ation

ike
Author Commented:
randDisp works fine, its at the bottom of the code snippet, simply returns floats in the range of the parameter passed (with a preference on positive values).
Commented:
I sort of got this working better....

The issue was the test for "doneness"....

this test
if ((xMidpnt - xMin) > 1) { // Subdivide if new midpoint length will be > 1

didn't return values in the expected range

it returned large intervals.... and it was never less than zero.

peace,
-john

Commented:
ie, xMin was always positive.. and so was xMidpiont....

the sides of the square were monotonically increasing... so in the interval 1 to 10 for instnace... xMid is 5, xMin is 1....

but midpoint never crosses 1.. so interval is big...

-john

Author Commented:
Jaggedness was being caused by overwriting edge midpoints that had already been computed.

Solution would be to code snippet to something like:

if(!yMap[xMin][zMidpnt])
yMap[xMin][zMidpnt] = (A + C + E) / 3 + randDisp(disp);  // F
if(!yMap[xMidpnt][zMin])
yMap[xMidpnt][zMin] = (A + B + E) / 3 + randDisp(disp);  // G
if(!yMap[xMax][zMidpnt])
yMap[xMax][zMidpnt] = (B + D + E) / 3 + randDisp(disp);  // H
if(!yMap[xMidpnt][zMax])
yMap[xMidpnt][zMax] = (C + D + E) / 3 + randDisp(disp);  // I

Checking to see if midpoints had been previously set.

Experts Exchange Solution brought to you by