Tempelmodelle in der Irrlicht Engine

by Jan Seyler und Nithi Rungtanapirom

Umsetzung

Teil 1: GIS-Daten

Diesen Teil der Umsetzung werde ich sehr kurz halten, da alles Nötige hierzu bereits in Kapitel "Hightmap aus Punktwolken -> Anwendung" beschrieben ist. Wir haben genau den dort beschriebenen Code als Anfang unseres Programms verwendet, wobei wir jedoch den Output als Tiff ausgeblendet haben.

Teil 2: Konfiguration

Zuerst noch eine Anmerkung:

Für unser Programm haben wir uns auf das Tempelgebiet beschränkt, dies kommt daher, dass der Rest des Gebiets, zu dem Daten vorliegen, keinerlei Erhöhungen oder ähnliches bietet und somit der wesentlich erhöhte Rechenaufwand nicht zweckmäßig ist.

Nun kommen wir zum grafischen Teil des Programms. Wir sind momentan auf folgendem Stand: Die erforderlichen GIS-Daten wurden eingelesen und alle fehlenden Höhen sind berechnet. D.h. das Array buffer16[] ist vollständig mit Höhen gefüllt. Nun müssen wir diese Höhen in die Engine einbringen.
Wir haben zuerst den Ansatz verfolgt den Zwischenschritt über eine Heightmap zu gehen, da Irrlicht in der Größe der einzulesenden Heightmap sehr beschränkt ist (das Maximum, das vollständig funktioniert hat, was 300x300 px) ist die Auflösung in Relation zur Terraingröße unzureichend.
Das ganze wird durch folgenden Screenshot verdeutlicht:
Problem Heightmap_Auflösung
Man erkennt, dass in dieser Skalierung die Höhen nicht wirklich sichtbar sind und das Gelände auch nicht richtig skaliert ist. Das Problem ergibt sich daraus, dass eine genaue Heightmap für das Gebiet eine Größe von 24000 x 21000 px bzw nur für das Tempelgebiet 8000 x 8000 px hat. Wenn man diese Heightmap nun zuerst auf 200 x 200 px runterrechnen muss und anschließend wieder auf eine einigermaßen realistische Dimension hochskalieren muss, sind weiterhin 70% der Daten verloren.
Für große Gelände ab ca. 1000 x 1000 px empfiehlt sich immer die Lösung, die wir gewählt haben, zu der ich gleich kommen werde.
Zunächst noch ein weiteres Defizit von Irrlicht, das man auf obigen Screenshot gut erkennen kann: Das LOD ist wesentlich zu grob (man betrachte den Berg links im Bild)! Soweit ich weiß, arbeiten die Entwickler bereits an einem dynamischen LOD, was auch sehr wünschenswert ist. Das ist ein weiteres Problem, warum es nicht sinnvoll ist, das ganze Terrain auf einmal zu rendern bzw. die Sichtweite zu groß zu wählen. Niemand möchte Berge und andere Details, die von einem zum nächsten Moment ihre Form komplett von einer Ecke zu 5 Ecken ändern. Um dieses Problem zu lösen, haben wir uns für die Lösung "Nebel" entschieden, was wiederrum andere Probleme aufwirft, dazu später mehr.
Nun zur Visualisierung großen Terrains.

Wir haben zunächst einmal ein Startfenster konstruiert, welches den Benutzer diverse Auswahlmöglichkeiten erlaubt.

MyEventReceiver receiver;
IrrlichtDevice *device;
video::E_DRIVER_TYPE drivertype = video::EDT_OPENGL;
bool fulscreen = false;

// setup

device = createDevice(drivertype, core::dimension2d<s32>(500, 300), 32, fulscreen, false, false, 0);

video::IVideoDriver* driver = device->getVideoDriver();
gui::IGUIEnvironment* guienv = device->getGUIEnvironment();

ITimer* timer = device->getTimer();

gui::IGUIFont* font = guienv->getFont("media/cityblueprint.bmp");
guienv->getSkin()->setFont(font);

Dies ist die Initialisierung für Irrlicht. Es wird u.a. der Grafiktreiber und die Fenstergröße festgelegt. Desweiteren definieren wir, wie in jedem Irrlicht-Programm zuerst einmal driver, guienv und timer. Zusätzlich übergeben wir der GUI noch eine eigene Schriftart.

device->setWindowCaption(L"Angkor Project");

guienv->addStaticText(L"Cambodia", core::rect<s32>(20,20,250,40), false);
guienv->addStaticText(L"created: 28.07.2008", core::rect<s32>(20,40,250,60), false);
guienv->addStaticText(L"author: Nithi und Jan", core::rect<s32>(20,60,250,80), false);

gui::IGUIButton *drtbutton = guienv->addButton(core::rect<s32>(170,70,330,100), 0, 101, L"OpenGl");
gui::IGUIButton *fscrbutton = guienv->addButton(core::rect<s32>(170,110,330,140), 0, 101, L"WINDOW");
gui::IGUIButton *contbutton = guienv->addButton(core::rect<s32>(170,150,330,180), 0, 101, L"CONTINUE");
gui::IGUIButton *extbutton = guienv->addButton(core::rect<s32>(170,190,330,220), 0, 101, L"EXIT");

bool exit = true;

drivertype = video::EDT_OPENGL;

u32 pause = 0;

Das eben erwähnte Startfenster wird hier nun eingerichtet. Zuerst setzen wir den Titel des Fensters. Danach haben wir etwas Text eingefügt, so dass man weiß, welches Objekt gestartet ist.
addStaticText() nimmt eigentlich sieben Argumente, wir ändern nur drei und die restlichen vier bleiben somit auf dem Standard.

Das erste Argument ist der Text, der angezeigt werden soll. Man kann dies später mit SetText() wieder abändern.
Das zweite Argument ist die Position des Textes, welche mit Hilfe eines Rechtecks übergeben wird. Um das Rechteck zu setzten, gibt man die obere linke Ecke und die rechte untere an. Man muss nur beachten, dass der Ursprung des zweidimensionalen Koordinatensystems in der oberen linken Ecke des Fensters liegt.
Das dritte Argument ist border, den man true oder false setzen kann. Hiermit entscheidet man, ob der Text einen 3D Rahmen hat.
Die letzten vier Argumente sind: wordWrap (Zeilenumbruch Standard="true"), parent (Standard=0), id (Standard= -1), fillBackground (Standard="false")

Anschließend fügen wir die Buttons für die Auswahlmöglichkeiten hinzu.
addButton() nimmt fünf Argumente.
Das erste Argument rectangel bestimmt wie eben die Position und Größe des Buttons
Das zweite und dritte Argument sind parent und id, wobei parent auf NULL gesetzt wird. Allerdings übergeben wir eine id, um später besser darauf zugreifen zu können.
Das vierte Argument übergibt den Text, der auf dem Button stehen soll.
Das letzte Argument, welches wir nicht belegen, setzt den Text, der im Tooltip gezeigt wird, aber haben keine Tooltips gesetzt, weil die Buttons leicht zu verstehen sind.

 while(device->run())
{
if(timer->getRealTime() > pause) pause = 0;

if(!pause)
{
if(drtbutton->isPressed())
{
pause = timer->getRealTime() + 300;

if(drivertype == video::EDT_OPENGL)
{
drivertype = video::EDT_DIRECT3D9;
drtbutton->setText(L"DirectX 9");
}
else
{
drivertype = video::EDT_OPENGL;
drtbutton->setText(L"OpenGl");
}
}

if(fscrbutton->isPressed())
{
pause = timer->getRealTime() + 500;

if(fulscreen)
{
fulscreen = false;
fscrbutton->setText(L"WINDOW");
}
else
{
fulscreen = true;
fscrbutton->setText(L"FULLSCREEN");
}

}
}

if(contbutton->isPressed())
{
exit = false;
device->closeDevice();
}

if(extbutton->isPressed())
{
exit = true;
device->closeDevice();
}

driver->beginScene(true, true, video::SColor(255,100,120,100));
guienv->drawAll();
driver->endScene();
}
device->drop();

if(exit) return 0;

In diesem Teil legen wir fest, was geschieht, wenn man die Buttons drückt. Das ganz geschieht natürlich während das device läuft, also während das Fenster angezeigt wird. Die Abfrage, ob ein button gedrückt ist oder nicht, wird wesentlich erleichtert dank der Irrlicht-Klasse irr::gui::IGUIButton. Hierbei kann man isPressed() abfragen.
Wenn man den contbutton drückt wird das nächste Fenster aufgerufen.
Nun beginnt die eigentliche Ausführung der Terrainvisualisierung.

 device = createDevice(drivertype, core::dimension2d<s32>(800, 600), 32, fulscreen, false, false, &receiver);
device->setWindowCaption(L"Angkor Project");

driver = device->getVideoDriver();
scene::ISceneManager* smgr = device->getSceneManager();
guienv = device->getGUIEnvironment();

font = guienv->getFont("media/cityblueprint.bmp");
guienv->getSkin()->setFont(font);

timer = device->getTimer();

// set attributes of terrain
s32 terrainWidth = (s32)(width);
s32 terrainHeight = (s32)(length);
s32 meshSize = 100;
f32 tileSize = 1;
core::vector3df terrainPos(0.0f, 0.0f, 0.0f);

// setup camera
MCameraFPS camera(smgr);
camera.setNearValue(0.1f);
camera.setFarValue(tileSize*meshSize/2);
camera.setPosition(core::vector3df(terrainWidth*tileSize/2, 0, terrainHeight*tileSize/2) + terrainPos);

// fog
driver->setFog(video::SColor(255,100,101,140), true, tileSize*meshSize/4, tileSize*(meshSize-4)/2, 0.05f);

Die Rückgaben aus den Buttons, werden nun für das zweite Device übergeben. Und das Fenster wird, genau wie das auch beim ersten geschehen ist, initialisiert. Um nicht die alten Werte und Einstellungen zu übernehmen, werden die Zeiger nocheinmal neu gesetzt.
Daran anschließend wird das Terrain eingerichtet. Die Höhe und Breite übernehmen wir so, wie wir es beim Einlesen der GIS-Daten gesetzt haben. Weil wir das auszulesende Gebiet verkleinert haben, gilt:

width=8000/r, height=8000/r

meshSize und tileSize sind die Eigenschaften, die für das Tiled Terrain wichtig sind. Was die Argumente beschreiben, zeigen die Variablennamen. terrainPos setzt den Nullpunkt des zukünftigen Terrains.
Danach initialisieren wir noch unsere Kamera, um diese zu verstehen, empfehle ich MCamera.h anzusehen, wie es im Source vorhanden ist. Im Grunde ist es nur die Irrlicht eigene Kamera, nur mit einigen Modifikationen, die z.B. ein 360° Blickfeld auch nach Oben und Unten ermöglicht. Zuletzt habe ich noch den Nebel gesetzt, damit man das Berechnen der Tiles beschränken, indem man nicht alles auf einmal Rendern muss, sondern nur 10 Tiles um die Position der Kamera. Somit werden deutlich Ressourcen gespart und die FPS ist wesentlich höher. Ein weiterer netter Nebeneffekt ist, dass das LOD wegfällt. Also kein nerviges und unschönes Popping des Geländes mehr.

Durch den Nebel ergibt sich allerdings ein neues Problem. Irrlicht kann Nebel und Skyboxen nicht vereinen und zeigt Nebel nur auf den Gegenständen, nicht aber vor den Skyboxen an. Dies ist auch einfach zu erklären: OpenGL rechnet Nebel abhängig von der Entfernung des Gegenstandes einfach als eine prozedurale Vermischung der Nebelfarbe mit den Farben des Objekts. Die Skybox lieg nun aber so weit entfernt, dass man entweder nur Nebel sehen könnte, wodurch man sich die Skybox natürlich sparen kann, oder eben die Skybox ohne Nebel. OpenGL entscheidet sich für die zweite Variante. Eine Lösung dieses Problems werde ich in "Ausblick" vorschlagen.

// light
scene::ILightSceneNode* light = smgr->addLightSceneNode(0, core::vector3df(0,0,0),
video::SColorf(255, 255, 255, 255), 1000.0f);
video::SLight ldata = light->getLightData();
ldata.AmbientColor = video::SColorf(0.2f,0.2f,0.2f);
ldata.DiffuseColor = video::SColorf(1.f,1.f,1.f);
ldata.Type = video::ELT_DIRECTIONAL;
ldata.Position = core::vector3df(-10.f,5.f,-5.f);
light->setLightData(ldata);

// create test scene node
scene::IAnimatedMesh* testmesh = smgr->getMesh("media/test.3ds");
scene::IAnimatedMeshSceneNode* test = smgr->addAnimatedMeshSceneNode( testmesh );
test->setMaterialFlag(video::EMF_LIGHTING, false);
test->setMaterialFlag(video::EMF_FOG_ENABLE, true);
cout << "create test scene node complete! \n";

In diesem Teil habe ich das Licht initialisiert, denn ohne Licht wäre das gesamte Bild später schwarz. Dann haben wir noch eine Test Scene Node eingefügt, einfach, damit man sieht, ob es geht. Und sonst kommen wir nun zum Terrain:

Teil 3: Tiled Terrain

Für unser Terrain haben wir uns der Node von arras bedient, aber einiges geändert.
Da für die Node selbst so gut wie keine Dokumentation besteht, werde ich auf die wichtigsten Aspekte genau eingehen.

ShTlTerrainSceneNode* terrain = new ShTlTerrainSceneNode(smgr, terrainWidth, terrainHeight, tileSize, meshSize); terrain->drop();

Nun schreibe ich also die eigentliche Terrain Node. Wie man sieht, ist terrain nun Teil der Klasse ShTlTerrainSceneNode. Als Argumente übergibt man die zuvor gesetzten Attribute des Terrains.

Um die Node im einzelnen zu verstehen, empfiehlt es sich dringend ShTlTerrainSceneNode.cpp durchzuarbeiten.

terrain->setMaterialFlag(video::EMF_LIGHTING, true);
terrain->setMaterialFlag(video::EMF_FOG_ENABLE, true);

terrain->setPosition(terrainPos);
cout << "setup terrain scene node complete! \n";

// set terrain to render around camera
terrain->follow(camera.getNode());

Hier werden, wie bei einem normalen Terrain, die Eigenschaften für den Treiber übergeben. Danach wird bestimmt, wo der Ursprung des Terrains liegen soll. Schließlich noch, wie bereits angedeutet, wird die Position der Kamera ausgelesen, so dass das Terrain um die aktuelle Position herum gerendert werden kann. Somit ist natürlich die maximale Größe des Terrains nicht mehr beschränkt. Da wir aber zuvor zu allen Punkten noch eine Höhe berechnen mussten, waren wir wegen der Laufzeit gezwungen das Gebiet zu beschränken.

makeHills(terrain, buffer16, terrainWidth, terrainHeight);
free(buffer16);
free(water);

Hier wird das Terrain errechnet. Um den Aufruf zu verstehen, zeige ich den entsprechenden Teil von SomeFunctions.h:

#include <irrlicht.h>
using namespace irr;

#include"ShTlTerrainSceneNode.h"
#include"random.h"

void makeHills(ShTlTerrainSceneNode* node, uint16* buffer16, s32 terrainWidth, s32 terrainHeight)
{

// make hills
for(s32 j=0; j<terrainHeight; j++)
{
for(s32 i=0; i<terrainWidth; i++)
{
int n=(int)(j*terrainWidth+i);
node->setHeight(i, j, (f32)buffer16[n]/128);
}
}
}

Wie man sieht steckt nicht viel hinter der Funktion. Nachdem die Höhenmatrix berechnet ist, muss man dem Programm nur noch zu jeder Koordniate die richtige Höhe via SetHeight() übergeben. Hier gehen wir, genau wie beim Beschreiben von buffer16[] wieder zeilenweise vor.

terrain->smoothNormals();

In diesem Schritt wird das Gelände ein wenig abgerundet.

randomizeColor(terrain);
randomizeUV(terrain);

Diese beiden Schritte lassen das Terrain wesentlich realistischer aussehen. randomizeColor() gibt dem Terrain ein zufällige Farbe, abhängig von der Höhe und randomizeUV() erstellt eine zufällige Anordung der UV-Vektoren, was gleichzeitig Texturen zufällig verteilt. Man kann dem Terrain beliebig viele Texturen übergeben und diese werden dann zufällig verteilt.

terrain->setMaterialTexture(1, driver->getTexture("media/dirt_x_43.bmp"));

driver->beginScene(true, true, video::SColor(255,100,101,140));
guienv->drawAll();
driver->endScene();

// movement speed per second
f32 speed = 10.0f * tileSize;

// time
f32 time = 0;
f32 oldtime;
// time between frames
f32 deltaT;

wchar_t tmp[1024];

bool help = false;

guienv->clear();

Als letztes wird nur noch die Detailmap, die Bewegungsgeschwindigkeit geschrieben. Außerdem brauchen wir, um später die FPS angeben zu können noch ein paar Variablen.

Teil 4: Wasserflächen

Um das Terrain realistischer wirken zu lassen, haben wir überall dort, wo in Realität auch Wasserflächen sind, Wasser gesetzt. Im Code sieht das ganze wie folgt aus:

//add animated water

scene::IAnimatedMesh* mesh1 = smgr->addHillPlaneMesh("water1#",
core::dimension2d<f32>(50,50),
core::dimension2d<u32>(20,10), 0, 0, //Größe (1000,500)
core::dimension2d<f32>(0,0),
core::dimension2d<f32>(20,20));

scene::ISceneNode* node1 = smgr->addWaterSurfaceSceneNode(mesh1->getMesh(0), 0.25f, 300.0f, 0.75f);
node1->setPosition(core::vector3df(1085,9.25f,1040)); //Mittelpunkt
node1->setMaterialTexture(0, driver->getTexture("Media/water_refl.jpg"));
node1->setMaterialTexture(1, driver->getTexture("Media/water.jpg"));

node1->setMaterialType(video::EMT_REFLECTION_2_LAYER);
node1->setMaterialFlag(video::EMF_FOG_ENABLE, true);

Wir haben zwar mehrere water meshes hinzugefügt, da aber alle genau gleich funktionieren, erkläre ich nur das erste.
Zuerst erschaffen wir ein leeres, flaches, animiertes Mesh, dem man dann noch Eigenschaften zuordnet. Zuallererst den Namen, den das Mesh später trägt. Die nächsten beiden Argumente ergeben zusammen die Größe des Meshs. Zuerst legt man die tileSize in x- und y-Richtung fest. Danach gibt man an, wie oft die Tiles in x- und y-Richtung jeweils wiederholt werden soll, also den tileCount. Insgesamt kann man dann die Größe des Meshs in jede Richtung ausrechnen, indem man tileSize * tileCount rechnet.
Nun könnte man noch das Material und die Hügelhöhe ändern, aber wir wollen ja ein leeres, flaches Mesh und lassen somit diese Werte auf NULL.
Das nächste Argument gibt die Anzahl an Hügeln an, die natürlich NULL sein soll. Schließlich legt man noch fest, wie oft die Textur in jede Richtung wiederholt werden soll, dies ist später der entscheidende Faktor für die Ressourcenaufwändigkeit des Wassers.
Nachdem wir nun ein leeres Mesh haben, wollen wir daraus natürlich noch Wasser machen. Dies tun wir in den nächsten Schritten. addWaterSurfaceSceneNode hat hierfür acht Argumente, von denen wir allerdings nur vier vom Standartwert her ändern.
Das Erste besagt, welches Mesh zu Wasser animiert werden soll. Wir wählen hier natürlich unser soeben erstelltes mesh1. Zweitens sagen wir, wie hoch die Wellen sein sollen, da wir kein Meer oder Ähnliches erstellen wollen, haben wir diese sehr niedrig gewählt. Nun folgen zwei weitere Argumente für die Wellen. Zuerst die Geschwindigkeit (niedrigere Zahl = höhere Geschwindigkeit) und dann die Länge der Wellen. Die letzten Argumente wären parent, id, position, rotation und scale diese ändern wir aber nicht.
Damit das Wasser schließlich an der gewollten Stelle ist, übergibt man nun noch den Mittelpunkt des Meshs.
Nun bearbeitet man noch das Material, übergibt also die Texturen und die Eigenschaften.

Teil 5: Radar

Da wir ein so großes Terrain haben, haben wir uns entschlossen, zur besseren Orientierung ein Radar einzubauen, Dies ist auch nicht weiter kompliziert, war aber schwerer als gedacht zu implementieren. Zuerst wollte ich alles in die While(device->run) Schleife einbauen, was aber die FPS von 180 auf 6 gedrückt hat. Also mussten wir uns eine andere Lösung überlegen und die sieht wie folgt aus:

 video::ITexture* dot=driver->getTexture("Media/dot.tga");
video::ITexture* map=driver->getTexture("Media/Map.jpg");

// Add Terrainradar

guienv->addImage(map,core::position2d<s32>(590,10));
//Map for radar
gui::IGUIElement* point=guienv->addImage(dot,core::position2d<s32>(590+100-2,10+100-2));
point;


while(device->run())
{...

Hier wird zuerst einmal die Karte und der Punkt in die Mitte der gesetzten Karte in das Fenster gelegt. Dabei nimmt addImage() lediglich zwei Argumente und zwar die Textur des geladenen Bildes und die Position der oberen linken Ecke des Bildes also 2D-Vector. Schließlich muss der Punkt sich nur noch genauso bewegen, wie es die Kamera tut:

point->setRelativePosition(core::rect<s32>(((s32)camera.getPosition().X/10)+590-2,210-((s32)camera.getPosition().Z/10)-2,((s32)camera.getPosition().X/10)+590+2,210-((s32)camera.getPosition().Z/10)+2));
}

Hierbei ist die eigene Kameraklasse wieder von Vorteil, da man die Position sehr einfach abfragen kann.
setRelativePosition() erwartet nur ein einziges Argument und zwar die neue Position als Rechteck. Dies geschieht genauso wie immer, wenn man ein Rechteck angibt (obere linke Ecke und untere rechte).

Teil 6: Waternode und Scanline-Algorithmus

Zuletzt biete ich einen Nachtrag, da wir diesen Teil erst am Ende erstellt haben. Wie oben erwähnt, wollten wir überall, wo in der Realität Wasser ist, auch in unserer Simulation Wasser setzen. Hierfür bietet es sich an, noch einmal die Satellitendaten zu verwenden. Das haben wir also auch gemacht:
In diesem Teil wird die Shapefile "Waterbody_Poly" gelesen, um die Lage des Wassergeländes zu bestimmen. Wir beginnen erstmal damit, die Shapedatei zu öffnen:

/* ################### Part 5: A ################### */
/* reading shapefiles 'Waterbody_Poly' */

cout << "reading Waterbody_Poly..." << endl;
startTime = clock();
shape=SHPOpen("WaterBody_Poly","rb");
SHPGetInfo(shape,&entities,&shptype,NULL,NULL);
cout << "Number of entities = " << entities << endl;
cout << "Shapetype: " << shptype << endl;

if( NULL == (water = (uint16*) malloc(pixelNumber*2)) )
{
fprintf(stderr, "Could not allocate enough memory for storing water!\n");
exit(EXIT_FAILURE);
}
for(i=0;i<entities;i++) water[i]=0;

Hierzu wurde das Array "water" eingerichtet. Jedes Element wird den Wert 0 oder 0xffff speichern, je nachdem, ob sich an derjenigen Stelle Wasser befindet. Nun kommen wir zu dem Scanline-Algorithmus:

 Dieser Algorithmus wird normalerweise verwendet, um ein Polygon zu füllen. Weitere Erklärungen findet sich z.B. auf Wikipedia.de. Das Polygon kann man hier als Wassergelände ansehen, und die Frage, ob ein Punkt innerhalb des Polygons liegt (d.h. ob der gegebene Pixel gefärbt werden soll), entspricht der Frage, ob sich an dieser Stelle Wasser befindet.

 Dazu wenden wir diesen Algorithmus auf jedes einzelne Objekt an. Ein Objekt besteht aus einem oder mehreren Polygonen. Im letzteren Fall besagt nParts die Anzahl der Polygone. Das erste Polygon ist das größte Polygon des Objekts, die anderen sind Polygone, die innerhalb des ersten Polygons liegen und somit "Löcher" im Polygon bilden. Den Beginn des neuen Polygons innerhalb des Objekts erkennt man an dem Array object->panPartStart[]. Die x- und y-Koordinaten der Eckpunkte der Polygone werden in padfX[] bzw. padfY[] gespeichert. Jedes Polygon in diesem Objekt beginnt und endet mit demselben Eckpunkt.

/* ################### Part 5: B ################### */
/* scanline algorithm */
cout << "applying scanline algorithm to each polygon... ";
for(i=0;i<entities;i++)
{
SHPObject* object=SHPReadObject(shape,i);
j=0;
minx=fx(object->dfXMin,min[0]); maxx=fx(object->dfXMax,min[0]);
miny=fy(object->dfYMax,max[1]); maxy=fy(object->dfYMin,max[1]);
// check if the polygon is inside the observed area.
if(miny <= length && maxy >= 0 && minx <= width && maxx >= 0)

Da wir uns auf das Tempelgebiet einschränken, können wir alle Polygone außerhalb dieses Gebiets ausschließen, ansonsten wird der Scanline-Algorithmus verwendet.

{
if( NULL == (x = (double*) malloc((object->nVertices)*sizeof(double))))
{
fprintf(stderr, "Could not allocate enough memory for storing all vertices!\n");
system("PAUSE");
exit(EXIT_FAILURE);
}
if( NULL == (y = (double*) malloc((object->nVertices)*sizeof(double))))
{
fprintf(stderr, "Could not allocate enough memory for storing all vertices!\n");
system("PAUSE");
exit(EXIT_FAILURE);
}
for(count=0; count<object->nVertices; count++)
{
x[count]=fx(object->padfX[count],min[0]);
y[count]=fy(object->padfY[count],max[1]);
}
if((int)miny<0) beginy=0; else beginy=(int)miny;
if((int)minx<0) beginx=0; else beginx=(int)minx;

Die x- und y-Koordinaten aller Eckpunkte werden in dem Feld x[] bzw. y[] gespeichert. Der bearbeitete Bereich für jedes einzelne Polygon ist gerade der Durchschnitt vom Tempelbereich und vom Polygonbereich (daher die Variablen beginx und beginy, und die untere Schleifen-Abbruchbedingung).

for(zeile=beginy; zeile<maxy+1 && zeile<length; zeile++) // Schleife nach Y-Koord.
{
// Bestimme alle schnittpunkte der Zeile mit allen Kanten des Polygons.
if( NULL == (p = (double*) malloc((object->nVertices-object->nParts)*sizeof(double))))
{
fprintf(stderr, "Could not allocate enough memory for storing all vertices!\n");
system("PAUSE");
exit(EXIT_FAILURE);
}
j=0; count=0;
for(k=0;k<object->nVertices-1;k++)
{
if(j<object->nParts-1 && k==object->panPartStart[j+1]-1) {j++;}
else if((y[k]>=zeile && y[k+1]<=zeile) || (y[k]<=zeile && y[k+1]>=zeile))
{
if(y[k]==zeile)
{
bool before, after;
if(k!=object->panPartStart[j]) before = (y[k-1]<zeile);
else if(j==object->nParts-1) before = (y[object->nVertices-2]<zeile);
else before = (y[object->panPartStart[j+1]-2]<zeile);
after = (y[k+1]>zeile);
if(before ^ after){ p[count]=x[k]; count++; }
}
else
{
p[count] = (zeile-y[k])*(x[k+1]-x[k])/(y[k+1]-y[k])+x[k];
count++;
}
}
}
qsort(p,count,sizeof(double),lessthan);

 

Nun kommen wir zum Hauptschritt des Algorithmus. Zuerst sind die Schnittpunkte jeder Bildzeile mit den Kanten des Polygons zu bestimmen. Dies wird in das Feld p[] geschrieben. Von Mathematik her ist ein Schnittpunkt gerade der von der Geraden y=zeile und der Strecke mit den Eckpunkten (x[k],y[k]) und (x[k+1],y[k+1]). Dies bedeutet, liegt zeile nicht zwischen y[k] und y[k+1], so gibt es keinen Schnittpunkt der Bildzeile mit der Kante. Ansonsten berechnet sich die x-Koordinate des Schnittpunkts (y-Koordinate ist wegen y=zeile bekannt!) aus der Formel (p[count]-x[k])/(zeile-y[k])=(x[k+1]-x[k])/(y[k+1]-y[k]) (die Kehrwerte der beiden Steigungen müssen gleich sein).

Allerdings sind noch ein paar Spezialfälle zu behandeln, nämlich wo y[k]=zeile. Hier müssen wir genauer betrachten, ob sich die Farbe der betrachteten Zeile an dieser Stelle ändern soll. Wird beispielsweise dieser Schnittpunkt genau einmal gespeichert, so ändert sich die Farbe an diesem Punkt. Wird dieser Schnittpunkt aber genau zweimal gespeichert, so ändert sich die Farbe nicht(Genaueres wird später erklärt). Danach werden alle Schnittpunkte sortiert.

// Dann alles berechnen
for(spalte=beginx; spalte<width && spalte<maxx+1; spalte++) // Schleife nach X-Koord.
{
int left=0;
for(j=0; j<count; j++){ if(p[j]<=spalte) left++; else break; }
if(left%2 == 1) water[zeile*width+spalte]=0xffff;
}
free(p);
}
free(x); free(y);
SHPDestroyObject(object);
}
}
cout << "done\n";
SHPClose(shape);

Jetzt kommen wir zum Schluss des Scanline-Algorithmus. Für jeden Punkt auf der Bildzeile wird abgezählt, wieviele Schnittpunkte mit allen Kanten des Polygons links von ihm liegt. Ist die Anzahl solcher Punkte ungerade, so wird dieser Punkt auf der Bildzeile gefärbt, d.h. er bekommt den Wert 0xffff zugewiesen. Rein theoretisch könnte man hier auch den anderen Wert wie z.B. 1 zuweisen, aber der Wert 0xffff entspricht der Farbe weiß, wenn wir später den "Watermap" erzeugen wollen.

/* ################### Part 5: C ################### */
/* optional: generating output 'watermap.tif' */

// Open file 'watermap.tif' for output
/* if((output = TIFFOpen("watermap.tif", "w")) == NULL)
{
fprintf(stderr, "Could not open outgoing watermap.tif\n");
exit(EXIT_FAILURE);
}

printf("writing watermap.tif ... ");

// Write the tiff tags to the file
TIFFSetField(output, TIFFTAG_IMAGEWIDTH, width);
TIFFSetField(output, TIFFTAG_IMAGELENGTH, length);
TIFFSetField(output, TIFFTAG_COMPRESSION, compression);
TIFFSetField(output, TIFFTAG_PLANARCONFIG, 1);
TIFFSetField(output, TIFFTAG_PHOTOMETRIC, PHOTOMETRIC_MINISWHITE);
TIFFSetField(output, TIFFTAG_BITSPERSAMPLE, 16);
TIFFSetField(output, TIFFTAG_SAMPLESPERPIXEL, 1);

// short and easy: write all the data into 'heightmap.tif'
if(TIFFWriteEncodedStrip(output, 0, water, pixelNumber*2 ) == 0)
{
fprintf(stderr, "Could not write image\n");
exit(42);
}

TIFFClose(output);
cout << "done\n";*/
endTime = clock();
printf("reading Waterbody_Poly and applying scanline algorithm... done: %.2f seconds\n", (float)(endTime-startTime) / CLOCKS_PER_SEC);
prgEnd = clock();
printf("\ntotal runtime for reading both shapefiles: %.2f seconds\n", (float)(prgEnd-prgStart) / CLOCKS_PER_SEC);
system("PAUSE");

In diesem Teil wird die entsprechende TIFF-Datei beschrieben, was aber ausgeblendet wurde. Schließlich wird der Anwender informiert, wie lange es gedauert hat, Waterbody_Poly auszulesen und den Scanline-Algorithmus darauf anzuwenden.