Générateur de planètes

Backup · Message par **Backup** » mar. 02/juil./2013 14:22

oui , faudra passer en Pb 5.20
car deja la fonction SetMeshData() n'existe pas en 5.11

GallyHC · Message par **GallyHC** » mar. 02/juil./2013 19:38

SPH > tu es en mode forage de vieux topic lollll

GallyHC

SPH · Message par **SPH** » mar. 02/juil./2013 19:43

GallyHC a écrit :SPH > tu es en mode forage de vieux topic lollll

GallyHC

C'est GG qui l'a déterré. Mais ca m'interesserais d'avoir la version 5.20b5

GallyHC · Message par **GallyHC** » mar. 02/juil./2013 20:21

Désole SPH

et pareil que mon dernier message pour GG alors.

GG · Message par GG » mer. 03/juil./2013 9:02

Salut JY,

J'ai déterré celui ci car je réfléchis à la conception d'un simulateur de système solaire...

GallyHC · Message par **GallyHC** » mer. 03/juil./2013 9:36

Gros projet en perspective

et c'été une boutade le mode forage

kelebrindae · Message par **kelebrindae** » jeu. 04/juil./2013 22:19

Ok, voici la version 5.20!
Les touches:
- curseur pour faire tourner la planète;
- T / R : augmenter / réduire le niveau de détail du mesh;
- W : afficher / cacher le maillage du mesh;
- La molette de la souris permet de zoomer.

N'oubliez pas de désactiver le debugger si ça rame trop (mais j'ai l'impression que c'est plus rapide qu'en 4.60; serait-ce dû à une optimisation de PB?)

Code : Tout sélectionner

; Author: Kelebrindae
; Date: december 27 2011
; PB version: v5.20

; -------------------------------------------------------------------------------------
; Purpose: Generate a random telluric planet
; -------------------------------------------------------------------------------------
; - To make make mapping easier and avoid distorsion to the poles, the mesh is actually
;   made of 6 planes that form a cube => all we need is 6 square heightmaps for the
;   relief, and 6 square textures.
; -------------------------------------------------------------------------------------
; Known bugs:
; - Visible seams between heightmaps
; - Visible seams between submeshes
; -------------------------------------------------------------------------------------

;- ---------- Controls ----------
;- Arrow keys to rotate
;- W : wireframe display
;- R/T : decrease/increase tesselation
;- Mouse wheel: zoom
;- ------------------------------

;- Constants
#PLANETSEED=987654321 ; change this to create a different planet
#WATERLEVEL=127
#RELIEFEXAGERATION=0.2 ; Change this to increase/decrease relief height
#USEAPIDRAWING = #False
#GRIDSIZE=256 ; heightmaps and textures size
#USEVORONOI = #True ; slower, but nicer

; Sides of the cube
Enumeration
  #TOPSIDE
  #BOTTOMSIDE
  #FRONTSIDE
  #BACKSIDE
  #LEFTSIDE
  #RIGHTSIDE
EndEnumeration

; Side of an image (to blend seams)
#HEIGHTMAP_TOP = 1
#HEIGHTMAP_BOTTOM = 2
#HEIGHTMAP_LEFT = 4
#HEIGHTMAP_RIGHT = 8


; Texture type
Enumeration
  #MAPTYPE_TERRAIN
  #MAPTYPE_ALTITUDE
  #MAPTYPE_SLOPE
  #MAPTYPE_TEMPERATURE
EndEnumeration


;- Data structures
Structure Vector3_struct ; used in vector maths (normalization)
  x.f
  y.f
  z.f
EndStructure

Structure VoronoiPoint2D_struct ; used to create the Voronoi Diagram (tectonic plates)
  x.i
  y.i
  dist.f
EndStructure

Structure geostat_struct ; used to store infos about terrains
  latitude.f
  slope.i
  temperature.f
  height.i
  probaterrain.f[10]
  isWater.b
  colour.i
EndStructure

Structure terrain_struct ; used to store specifications for each terrain type
  name.s
  heightAv.f
  heightVt.f
  slopeAv.f
  slopeVt.f
  temperatureAv.f
  temperatureVt.f
  colour.i
EndStructure

Structure vertex_struct
  x.f
  y.f
  z.f
  nx.f
  ny.f
  nz.f
  colour.i
  u.f
  v.f
EndStructure

Structure triangle_struct
  v1.i
  v2.i
  v3.i
 
  *ptrV1.vertex_struct
  *ptrV2.vertex_struct
  *ptrV3.vertex_struct
EndStructure

; API drawing
CompilerIf #USEAPIDRAWING = #True
  Structure myBITMAPINFO
    bmiHeader.BITMAPINFOHEADER
    bmiColors.RGBQUAD[1]
  EndStructure
CompilerEndIf


;- Global definitions
Global CameraMode.b
Global anglex.f = 0,angley.f = 0


Global Dim diamond.b(#GRIDSIZE,#GRIDSIZE)
Global Dim voronoi.b(#GRIDSIZE,#GRIDSIZE)

; There are 6 heightmaps and geostats arrays: one for each side of the cube
Global Dim heightmap.i(6,#GRIDSIZE,#GRIDSIZE)
Global Dim newhmap.i(6,#GRIDSIZE,#GRIDSIZE)
Global Dim geostat.geostat_struct(6,#GRIDSIZE,#GRIDSIZE)

;- Terrain types
Global Dim terrain.terrain_struct(10)

; Temperatures for poles and equator (Celsius)
Global tpol.f = -15,teq.f = 50,name.s

; Vertices which are on the heightmaps' borders
Structure borderVert_struct
  *ptrVert.vertex_struct
  u.f
  v.f
EndStructure

; Mesh info for each side of the cube
Structure planetSide_struct
  numTexture.i
  numMaterial.i
 
  List topBorder.borderVert_struct()
  List bottomBorder.borderVert_struct()
  List leftBorder.borderVert_struct()
  List rightBorder.borderVert_struct()
 
  List vertex.vertex_struct()
  List triangle.triangle_struct()
EndStructure
Global Dim planetSide.planetSide_struct(6)

; Planet entity
Global planetEntity.i, tesselation.i = 5


EnableExplicit

;************************************************************************************
;-                                 ---- Macros ----
;************************************************************************************

Global SQR_sqrhalf.f,*SQR_sqrptr,SQR_sqrInt.i
; Fast inverse square root, Romero style
Macro FAST_INVSQR(x,result)
  result = x
  SQR_sqrhalf = 0.5 * result
  *SQR_sqrptr = @result
  SQR_sqrInt  = PeekI(*SQR_sqrptr)
 
  ; The magic number is different in 64 bits
  CompilerIf #PB_Compiler_Processor = #PB_Processor_x64
    PokeI(*SQR_sqrptr,$5FE6EB50C7AA19F9 - ( SQR_sqrInt >> 1))
  CompilerElse ; all the others are 32 bits ?
    PokeI(*SQR_sqrptr,$5F375A86 - ( SQR_sqrInt >> 1)) ; or $5f3759df
  CompilerEndIf
 
  result=result*(1.5 - SQR_sqrhalf * result * result)
  result=result*(1.5 - SQR_sqrhalf * result * result) ; twice for greater precision
  result=result*(1.5 - SQR_sqrhalf * result * result) ; thrice is even better (but you can ditch it to squeeze a few extra cycles from your code)
EndMacro

; Vector3 normalization
Global VEC3_length.f
Macro VEC3_NORMALIZE(vector)
  FAST_INVSQR((vector\x*vector\x) + (vector\y*vector\y)+(vector\z*vector\z),VEC3_length)
  vector\x * VEC3_length
  vector\y * VEC3_length
  vector\z * VEC3_length
EndMacro

Macro VEC3_NORMALIZETOLENGTH(vector,magnitude)
  FAST_INVSQR((vector\x*vector\x) + (vector\y*vector\y)+(vector\z*vector\z),VEC3_length)
  VEC3_length*(magnitude)
  vector\x * VEC3_length
  vector\y * VEC3_length
  vector\z * VEC3_length
EndMacro

; Commpute a 3D point between two others (for tesselation)
Macro MIDDLEPOINT(result,pt1,pt2)
  result\x = pt1\x + (pt2\x - pt1\x) / 2
  result\y = pt1\y + (pt2\y - pt1\y) / 2
  result\z = pt1\z + (pt2\z - pt1\z) / 2
 
  result\u = pt1\u + (pt2\u - pt1\u) / 2
  result\v = pt1\v + (pt2\v - pt1\v) / 2
EndMacro

; Add a point in the list, only if it doesn't already exists. If it exists, return the existing point
Macro CREATENEWVERTEX(refPt,numNewPt)
 
  ; This uses the concatenated X,Y,Z coords as a map key to check for duplicate vertices (yeah, I know, it's ugly)
  ; ("Strf(...,2)" assumes that vertex coords within the 0.01 threshold are identical)
  numNewPt = antidup(StrF(refPt\x,2)+StrF(refPt\y,2)+StrF(refPt\z,2))-1
  If numNewPt = -1
    AddElement(newVertex())
    newVertex()\x = refPt\x
    newVertex()\y = refPt\y
    newVertex()\z = refPt\z
    newVertex()\nx = refPt\nx
    newVertex()\ny = refPt\ny
    newVertex()\nz = refPt\nz
    newVertex()\u = refPt\u
    newVertex()\v = refPt\v
    newVertex()\colour = refPt\colour
   
    numNewPt = ListIndex(newVertex())
    AddMapElement(antidup(),StrF(refPt\x,2)+StrF(refPt\y,2)+StrF(refPt\z,2),#PB_Map_NoElementCheck)
    antidup() = numNewPt+1
  EndIf
   
EndMacro

; Weld the vertices between two heightmaps.
; Ex: weld the bottom of the "Top" heightmap to the top of the "Front" heightmap
Macro WELDBORDERVERTICES(side1,border1,coord1,side1Sort,side2,border2,coord2)
  SortStructuredList(planetSide(side1)\border1#Border(),side1Sort,OffsetOf(borderVert_struct\coord1),#PB_Float)
  SortStructuredList(planetSide(side2)\border2#Border(),#PB_Sort_Ascending,OffsetOf(borderVert_struct\coord2),#PB_Float)
 
  FirstElement(planetSide(side2)\border2#Border())
  ForEach planetSide(side1)\border1#Border()   
    planetSide(side2)\border2#Border()\ptrVert\x = planetSide(side1)\border1#Border()\ptrVert\x
    planetSide(side2)\border2#Border()\ptrVert\y = planetSide(side1)\border1#Border()\ptrVert\y
    planetSide(side2)\border2#Border()\ptrVert\z = planetSide(side1)\border1#Border()\ptrVert\z

    NextElement(planetSide(side2)\border2#Border())
  Next planetSide(side1)\border1#Border()

EndMacro

;************************************************************************************
;-                                 ---- Procedures ----
;************************************************************************************

;- --- Heightmaps generation ---
;************************************************************************************
; Name: makeDiamondSquare
; Purpose: Creates an array of fractal noise (looks like Perlin) that can be used as
;          a heightmap
; Parameters:
;   - size of the array (a square of size x size)
;   - dispersion
;   - seed: albeit the noise is random, a seed always produces the same result
; Return-value: none, but the result is stored in the "diamond" array
;************************************************************************************
Procedure makeDiamondSquare(gridSize.i,dispersion.f,seed.i)

  Protected gridstep.i,mid.i,n.i
  Protected i.i,i1.i,i2.i
  Protected j.i,j1.i,j2.i,js.i
  Protected u.f,average.f
  Protected dispStep.f,min.f=999999,max.f=-999999,ratio.f
  Protected Dim temp.f(gridSize,gridSize)
  Protected t.i = ElapsedMilliseconds()

  RandomSeed(seed)

  gridstep = gridSize
  dispStep = dispersion

  ; main loop
  While gridstep > 1
    mid=gridstep / 2

    ; Diamond step - calculates new diamond corners from squares
    i=mid
    i1 = i - mid
    i2 = i + mid
   
    While i < gridSize
      j = mid
      j1 = j - mid
      j2 = j + mid
     
      While j < gridSize
        ; Average of surrounding points
        average = (temp(i1,j1)+temp(i1,j2)+temp(i2,j1)+temp(i2,j2))/4
       
        ; calculate random values between -1 and 1
        u = Random(16384)/8192 - 1
       
        ; Diamond value
        temp(i,j) = average + (u * dispStep)
        j + gridstep
        j1 + gridstep
        j2 + gridstep
      Wend
     
      i + gridstep
      i1 + gridstep
      i2 + gridstep
    Wend

    ; square Step - calculates new square corners from diamonds
    i = 0
    i1 = i - mid
    i2 = i + mid
    js = 0
   
    While i < gridSize
      js = mid - js ; toggle start values of j loop
      j = js
      j1 = j - mid
      j2 = j + mid
     
      While j < gridSize ;+1
        average = 0
        If i1 < 0  ; check for need to wrap around i value
          average + temp(i2,j) + temp(i2,j)
        Else
          If i2 > gridSize
            average + temp(i1,j) + temp(i1,j)
          Else
            average + temp(i1,j) + temp(i2,j)
          EndIf
        EndIf
       
        If j1<0  ; check for need to wrap around j value
          average + temp(i,j2) + temp(i,j2)
        Else
          If j2 > gridSize
            average + temp(i,j1) + temp(i,j1)
          Else
            average + temp(i,j1) + temp(i,j2)
          EndIf
        EndIf     
        average/4
       
        ; calculate random value between -1 And 1
        u = Random(16384)/8192 - 1
       
        temp(i,j) = average + (u * dispStep)
        temp(gridSize,j) = temp(0,j) ; copy opposite edge
        j + gridstep
        j1 + gridstep
        j2 + gridstep
       
      Wend
      If j = gridSize
        temp(i,j) = temp(i,0) ; copy opposite edge
      EndIf
      i + mid
      i1 + mid
      i2 + mid
    Wend

    dispStep/2
    gridstep/2
  Wend

; Keep values in Byte range (between 0 and 255)
  For i=0 To gridSize-1
    For j=0 To gridSize-1
      ; Min /max values
      If temp(i,j)<min
        min=temp(i,j)
      Else
        If temp(i,j)>max
          max=temp(i,j)
        EndIf
      EndIf   
    Next j
  Next i
 
  ratio=255.0/(max-min)
  For i=0 To gridSize-1
    For j=0 To gridSize-1
      n=(temp(i,j)-min)*ratio
      If n<0
        diamond(i,j)=0
      Else
        If n>255
          diamond(i,j)=255
        Else
          diamond(i,j)=n
        EndIf
      EndIf
    Next j
  Next i
 
  Debug "Diamond square: "+Str(ElapsedMilliseconds()-t)+" ms."

EndProcedure

;************************************************************************************
; Name: makeVoronoi
; Purpose: Creates a modified Voronoi diagram that can be used to simulate tectonic
;          plates
; Parameters:
;   - size of the array (a square of size x size)
;   - number of regions
;   - number of points in each region (a point is the "center" of a Voronoi plate)
;   - seed: albeit the diagram is random, a seed always produces the same result
; Return-value: none, but the result is stored in the "voronoi" array
;************************************************************************************
Procedure makeVoronoi(gridSize.i,nbregion.i,nbpts.i,seed.i)

  Protected NewList ftpoint.VoronoiPoint2D_struct()
  Protected Dim coef.f(nbregion*nbregion*nbpts)
  Protected i.i,j.i,k.i,n.i
  Protected min.f=999999,max.f=-999999,ratio.f
  Protected Dim temp.f(gridSize,gridSize)
  Protected t.i = ElapsedMilliseconds()
 
  RandomSeed(seed)
  coef(0)=-1
  coef(1)=-2
  coef(2)=0.5
 
  ; In each region of the grid, place some points
  For i=0 To nbregion-1
    For j =0 To nbregion-1
      n=Random(nbpts)
      For k=1 To n
        AddElement(ftpoint())
        ftpoint()\x=i*(gridSize/nbregion)+Random(gridSize/nbregion)
        ftpoint()\y=j*(gridSize/nbregion)+Random(gridSize/nbregion)
        ;Debug StrF(ftpoint()\x)+","+StrF(ftpoint()\y)
      Next k
    Next j
  Next i
 
  ; For each cell of the grid, compute distance to each point, sort points according to distance.
  ; => value = coef * distance
  For i=0 To gridSize-1
    For j=0 To gridSize-1
   
      ForEach ftpoint()     
        ftpoint()\dist = Sqr((ftpoint()\x-i)*(ftpoint()\x-i)+(ftpoint()\y-j)*(ftpoint()\y-j))
      Next
      SortStructuredList(ftpoint(),#PB_Sort_Ascending,OffsetOf(VoronoiPoint2D_struct\dist),#PB_Float)
     
      k=0
      ForEach ftpoint()
        temp(i,j) + coef(k)*ftpoint()\dist
        k+1
      Next     
     
      ; Min /max values
      If temp(i,j)<min
        min=temp(i,j)
      Else
        If temp(i,j)>max
          max=temp(i,j)
        EndIf
      EndIf   
     
    Next j
  Next i

  ; Keep values in Byte range (between 0 and 255)
  ratio=255.0/(max-min)
  For i=0 To gridSize-1
    For j=0 To gridSize-1

      n=(temp(i,j)-min)*ratio
      If n<0
        voronoi(i,j)=0
      Else
        If n>255
          voronoi(i,j)=255
        Else
          voronoi(i,j)=n
        EndIf
      EndIf

    Next j
  Next i
 
  Debug "Voronoi: "+Str(ElapsedMilliseconds()-t)+" ms."
EndProcedure

;************************************************************************************
; Name: MakeHeightmap
; Purpose: Mix a diamond-square array and a Voronoi diagram to obtain the final heightmap
; Parameters:
;   - number of the heightmap
;   - size of the array (a square of size x size); Diamond, Voronoi, and Heightmap must
;     have the same size
;   - use the diamond-square array or not
;   - use the Voronoi diagram or not
; Return-value:  none, but the result is stored in the "heightmap" array
;************************************************************************************
Procedure makeHeightmap(numHmap.i,gridSize.i,Diamond.b,Voronoi.b)

  Protected i.i,j.i,k1.i,k2.i
  Protected min.f=999999,max.f=-999999,ratio.f,temp.f
  Protected Dim temp.f(gridSize,gridSize)
  Protected t.i = ElapsedMilliseconds()
 
  ; for each point of the array:
  For i=0 To gridSize-1
    For j=0 To gridSize-1
   
      ; both arrays will be used
      If Voronoi=#True And Diamond=#True
        k1=(diamond(i,j) & 255)
        k2=(voronoi(i,j) & 255)
        temp(i,j) = k1*k2
      Else
        ; only diamond-square
        If Diamond=#True
          k1=(diamond(i,j) & 255)
          temp(i,j) = k1
        Else
          ; only Voronoi
          If Voronoi=#True
            k2=(voronoi(i,j) & 255)
            temp(i,j) = k2
          EndIf
        EndIf
      EndIf   
   
      ; Store Min / Max values
      If temp(i,j)<min
        min=temp(i,j)
      Else
        If temp(i,j)>max
          max=temp(i,j)
        EndIf
      EndIf
   
    Next j
  Next i
 
  ; The heightmap values must be between 0 and 255
  ratio = 255.0/(max-min)
  For i=0 To gridSize-1
    For j=0 To gridSize-1
      temp = (temp(i,j)-min)*ratio
      ;Debug temp
      If temp<0
        heightmap(numHmap,i,j)=0
      Else
        If temp>255
          heightmap(numHmap,i,j)=255
        Else
          heightmap(numHmap,i,j)=temp
        EndIf
      EndIf
    Next j
  Next i
 
  Debug "Make heightmap: "+Str(ElapsedMilliseconds()-t)+" ms."

EndProcedure

;************************************************************************************
; Name: eraseHeightmapSeam
; Purpose: Blends the seam between two heightmaps
; Parameters:
;   - number of the first heightmap
;   - side of the first heightmap to blend (top, bottom, left, or right)
;   - number of the second heightmap
;   - side of the second heightmap to blend (top, bottom, left, or right)
;   - size of the array (a square of size x size); both Heightmaps must
;     have the same size
;   - percentage of that size that will be blended
; Return-value:  none
;************************************************************************************
Procedure eraseHeightmapSeam(numHmap1.i,side1.b,numHmap2.i,side2.b,gridSize.i,percent.b)

  Protected i.i,j.i,width.i,seamwidth.i,color.i
  Protected x1.i,y1.i,x2.i,y2.i
  Protected h1.i,h2.i,h.i
  Protected ratio1.f,ratio2.f
  Protected t.i = ElapsedMilliseconds()
 
  width = gridSize-1
  seamwidth = ( width*percent/100 )/2
   
  For i=0 To seamwidth
 
    ratio1 = (i/2)/seamwidth
    ratio2 = 1.0 - ratio1

    For j=0 To gridSize-1
   
      Select side1
        Case #HEIGHTMAP_TOP
          x1 = j
          y1 = seamwidth-i
        Case #HEIGHTMAP_BOTTOM
          x1 = j
          y1 = width-seamwidth+i
        Case #HEIGHTMAP_LEFT
          x1 = seamwidth-i
          y1 = j
        Case #HEIGHTMAP_RIGHT
          x1 = width-seamwidth+i
          y1 = j
      EndSelect
     
      Select side2
        Case #HEIGHTMAP_TOP
          x2 = j
          y2 = seamwidth-i
          If side1 = #HEIGHTMAP_TOP Or side1 = #HEIGHTMAP_RIGHT
            x2=gridSize-1-j
          EndIf         
        Case #HEIGHTMAP_BOTTOM
          x2 = j
          y2 = width-seamwidth+i
          If side1 = #HEIGHTMAP_BOTTOM Or side1 = #HEIGHTMAP_LEFT
            x2=gridSize-1-j
          EndIf
        Case #HEIGHTMAP_LEFT
          x2 = seamwidth-i
          y2 = j
          If side1 = #HEIGHTMAP_BOTTOM Or side1 = #HEIGHTMAP_LEFT
            x2=gridSize-1-j
          EndIf
        Case #HEIGHTMAP_RIGHT
          x2 = width-seamwidth+i
          y2 = j
          If side1 = #HEIGHTMAP_TOP Or side1 = #HEIGHTMAP_RIGHT
            x2=gridSize-1-j
          EndIf         
      EndSelect
   
      h1=heightmap(numHmap1,x1,y1)
      h2=heightmap(numHmap2,x2,y2)

      h=h1*ratio2 + h2*ratio1     
      heightmap(numHmap1,x1,y1)=h
     
      h=h2*ratio2 + h1*ratio1
      heightmap(numHmap2,x2,y2)=h     
    Next j
   
  Next i
 
  Debug "Erase seams: "+Str(ElapsedMilliseconds()-t)+" ms."

EndProcedure


;************************************************************************************
; Name: ComputeGeoStats
; Purpose: Compute height, slope, temperature, ice/snow, terrain type and color for
;          each cell of an heightmap
; Parameters:
;   - number of the heightmap
;   - size of the array (a square of size x size)
;   - sea level (between 0-255)
;   - Temperature at poles (usually cold)
;   - Temperature at equator (usually hot)
;   - seed: albeit the map is random, a seed always produces the same result
; Return-value:  none, but the results are stored in the "geostat" array
;************************************************************************************
Procedure computeGeoStats(numHmap.i,gridSize.i,waterLevel.i,tpol.f,teq.f,seed.i)

  Protected heightAboveWater.i
  Protected equivHeightToTemp.f = -0.4; temperature loss when you climb up 1 unit
  Protected i.i,j.i,k.i
  Protected halfGridSize.i = gridSize/2
  Protected x.i,y.i,newx.i,newy.i
  Protected p.i,h1.i,h2.i,maxp.i
  Protected red1.i,red2.i,green1.i,green2.i,blue1.i,blue2.i,red.f,green.f,blue.f
  Protected exp.f,cumul.f,icegradient.f,icecolor.i
  Protected vector3.Vector3_struct
  Protected t.i = ElapsedMilliseconds()

  RandomSeed(seed)
 
  gridSize-1
  For i=0 To gridSize
    For j = 0 To gridSize
   
      ; Latitude
      If numHmap=#TOPSIDE Or numHmap=#BOTTOMSIDE
        vector3\x = i - halfGridSize
        vector3\y = halfGridSize
        vector3\z = j - halfGridSize
   
        VEC3_NORMALIZE(vector3)
       
        vector3\x * halfGridSize
        vector3\y * halfGridSize
        vector3\z * halfGridSize

        geostat(numHmap,i,j)\latitude = Cos(vector3\y/gridSize * #PI)
      Else
        vector3\x = i - halfGridSize
        vector3\y = j - halfGridSize
        vector3\z = halfGridSize
   
        VEC3_NORMALIZE(vector3)
       
        vector3\x * halfGridSize
        vector3\y * halfGridSize
        vector3\z * halfGridSize

        geostat(numHmap,i,j)\latitude = Cos(vector3\y/gridSize * #PI)
      EndIf
   
      ; Altitude
      geostat(numHmap,i,j)\height = heightmap(numHmap,i,j)
     
      ; Slope
      h1= geostat(numHmap,i,j)\height
      geostat(numHmap,i,j)\slope=0
     
      If i>0
        h2=heightmap(numHmap,i-1,j)
        p=Abs( h1 - h2 )
        If p>geostat(numHmap,i,j)\slope
          geostat(numHmap,i,j)\slope = p
        EndIf
      EndIf

      If i<gridSize
        h2=heightmap(numHmap,i+1,j)
        p=Abs( h1 - h2 )
        If p>geostat(numHmap,i,j)\slope
          geostat(numHmap,i,j)\slope = p
        EndIf
      EndIf

      If j>0
        h2=heightmap(numHmap,i,j-1)
        p=Abs( h1 - h2 )
        If p>geostat(numHmap,i,j)\slope
          geostat(numHmap,i,j)\slope = p
        EndIf
      EndIf

      If j<gridSize
        h2=heightmap(numHmap,i,j+1)
        p=Abs( h1 - h2 )
        If p>geostat(numHmap,i,j)\slope
          geostat(numHmap,i,j)\slope = p
        EndIf
      EndIf
      If geostat(numHmap,i,j)\slope > maxp
        maxp = geostat(numHmap,i,j)\slope
      EndIf

      ; Temperature
      heightAboveWater=geostat(numHmap,i,j)\height - waterlevel
      If heightAboveWater <0
        heightAboveWater = 0
      EndIf
      geostat(numHmap,i,j)\temperature = tpol + geostat(numHmap,i,j)\latitude*(teq-tpol) + equivHeightToTemp*heightAboveWater
     
      If geostat(numHmap,i,j)\height<waterlevel
        ; This point is under sea level
        geostat(numHmap,i,j)\isWater=2
     
        ; underwater points => shades of blue
        red1=0:green1=3:blue1=92
        red2=0:green2=130:blue2=220
       
        red=red1+(red2-red1)*(geostat(numHmap,i,j)\height/waterlevel)
        green=green1+(green2-green1)*(geostat(numHmap,i,j)\height/waterlevel)
        blue=blue1+(blue2-blue1)*(geostat(numHmap,i,j)\height/waterlevel)
     
        ; water depth causes temperature variation
        geostat(numHmap,i,j)\temperature-(geostat(numHmap,i,j)\height-waterlevel)*0.075
      Else       
        ; This point is above sea level
        geostat(numHmap,i,j)\isWater=0
       
        ; Probability computation for each terrain type
        cumul=0:red=0:green=0:blue=0
        For k=1 To 5
          exp = Abs(terrain(k)\temperatureAv - geostat(numHmap,i,j)\temperature)/terrain(k)\temperatureVt
          exp + Abs(terrain(k)\slopeAv - geostat(numHmap,i,j)\slope)/terrain(k)\slopeVt
          exp + Abs(terrain(k)\heightAv*waterLevel - geostat(numHmap,i,j)\height)/terrain(k)\heightVt
          geostat(numHmap,i,j)\probaterrain[k] = Exp(-exp)
         
          cumul + geostat(numHmap,i,j) \ probaterrain[k]
        Next k
       
        ; Mix terrains colors according to probability
        For k=1 To 5
          geostat(numHmap,i,j)\probaterrain[k] / cumul
          blue+Red(terrain(k)\colour)*geostat(numHmap,i,j)\probaterrain[k]
          green+Green(terrain(k)\colour)*geostat(numHmap,i,j)\probaterrain[k]
          red+Blue(terrain(k)\colour)*geostat(numHmap,i,j)\probaterrain[k]
        Next k
               
      EndIf
 
      geostat(numHmap,i,j)\colour = RGB(red,green,blue)
     
    Next j
  Next i
 
  ; Generate a few rivers
  j=0
  For i=1 To Random(6)
    Repeat
      x=Random(gridSize)
      y=Random(gridSize)
      j+1
    Until geostat(numHmap,x,y)\isWater<>2 Or j=1000
   
    Repeat
      geostat(numHmap,x,y)\colour = RGB(0,130,220)
      geostat(numHmap,x,y)\isWater=1
      h1=999
      newx=0:newy=0
      If x>0
        If geostat(numHmap,x-1,y)\isWater <>1 And geostat(numHmap,x-1,y)\height < h1
          newx=x-1:newy=y
          h1=geostat(numHmap,newx,newy)\height
        EndIf
      EndIf
      If x<gridSize
        If geostat(numHmap,x+1,y)\isWater <>1 And geostat(numHmap,x+1,y)\height < h1
          newx=x+1:newy=y
          h1=geostat(numHmap,newx,newy)\height
        EndIf
      EndIf
      If y>0
        If geostat(numHmap,x,y-1)\isWater <>1 And geostat(numHmap,x,y-1)\height < h1
          newx=x:newy=y-1
          h1=geostat(numHmap,newx,newy)\height
        EndIf
      EndIf
      If y<gridSize
        If geostat(numHmap,x,y+1)\isWater <>1 And geostat(numHmap,x,y+1)\height < h1
          newx=x:newy=y+1
          h1=geostat(numHmap,newx,newy)\height
        EndIf
      EndIf
      x=newx:y=newy
     
    Until geostat(numHmap,x,y)\isWater=2 Or (newx=0 And newy=0)
   
  Next i

  ; Add ice/snow, slope marking
  For i=0 To gridSize
    For j=0 To gridSize
      red=Red(geostat(numHmap,i,j)\colour)
      green=Green(geostat(numHmap,i,j)\colour)
      blue=Blue(geostat(numHmap,i,j)\colour)
   
      ; Ice / snow
      If geostat(numHmap,i,j)\temperature < 0
        If geostat(numHmap,i,j)\temperature < 0
          icegradient= -geostat(numHmap,i,j)\temperature / 5
        EndIf
       
        If icegradient>1
          icegradient=1
        EndIf         
        icecolor=247+Random(8)
       
        red=red*(1-icegradient) + icecolor*icegradient
        green=green*(1-icegradient) + icecolor*icegradient
        blue=blue*(1-icegradient) + icecolor*icegradient       
      EndIf
     
      ; Slopes are darker
      If geostat(numHmap,i,j)\isWater<2   
        red-geostat(numHmap,i,j)\slope * 3
        green-geostat(numHmap,i,j)\slope * 3
        blue-geostat(numHmap,i,j)\slope * 3       
      EndIf     
       
      ; RGB must stay between 0-255
      If red<0
        red=0
      Else
        If red>255
          red=255
        EndIf
      EndIf
      If green<0
        green=0
      Else
        If green>255
          green=255
        EndIf
      EndIf
      If blue<0
        blue=0
      Else
        If blue>255
          blue=255
        EndIf
      EndIf
       
      geostat(numHmap,i,j)\colour = RGB(blue,green,red)

    Next j
  Next i
 
  Debug "Compute geostats: "+Str(ElapsedMilliseconds()-t)+" ms."

EndProcedure

;************************************************************************************
; Name: MakeTerrainImage
; Purpose: Draw a map from heightmap and geostats arrays
; Parameters:
;   - number of the heightmap / geostats array
;   - size of the array (a square of size x size)
;   - Type of map to draw: "terrain","altitude","slope", "temperature"
; Return-value:  generated image number
;************************************************************************************
Procedure.i makeTerrainImage(numHmap.i,gridSize.i,imageType.b)
 
  Protected numImage.i
  Protected i.i,j.i,color.i,icecolor.i
  Protected red.f,green.f,blue.f,icegradient.f
  Protected t.i =ElapsedMilliseconds()
 
  numImage = CreateImage(#PB_Any,gridSize,gridSize)
 
  CompilerIf #USEAPIDRAWING = #True
    Protected hDC.i,hBmp.i,*mem.i,picl_X.i,picl_Y.i,picl_D.i
    Protected bmi.myBITMAPINFO,*pixel.LONG
 
    hDC=StartDrawing(ImageOutput(numImage))
    hBmp = ImageID(numImage)
    picl_X = ImageWidth(numImage)
    picl_Y = ImageHeight(numImage)+1
    picl_D = ImageDepth(numImage) 
   
    *mem = AllocateMemory(picl_X*picl_Y*4)
    bmi.myBITMAPINFO
    bmi\bmiHeader\biSize   = SizeOf(BITMAPINFOHEADER)
    bmi\bmiHeader\biWidth  = picl_X
    bmi\bmiHeader\biHeight = picl_Y
    bmi\bmiHeader\biPlanes = 1
    bmi\bmiHeader\biBitCount = 32
    bmi\bmiHeader\biCompression = #BI_RGB
    GetDIBits_(hDC,hBmp,1,picl_Y,*mem,bmi,#DIB_RGB_COLORS)
 
    *pixel = *mem
  CompilerElse
    StartDrawing(ImageOutput(numImage))
  CompilerEndIf
 
  For j=gridSize-1 To 0 Step-1
    For i=0 To gridSize-1
   
      Select imageType
     
        Case #MAPTYPE_ALTITUDE
          red = geostat(numHmap,i,j)\height
          green = red:blue = Red ; (greyscale)

         
        Case #MAPTYPE_TEMPERATURE
          If i>1 And i<gridSize-1 And j>1 And j<gridSize-1 And geostat(numHmap,i,j)\isWater>0 And (geostat(numHmap,i-1,j)\isWater=0 Or geostat(numHmap,i+1,j)\isWater=0 Or geostat(numHmap,i,j-1)\isWater=0 Or geostat(numHmap,i,j+1)\isWater=0)
            color=0
          Else
            red = (geostat(numHmap,i,j)\temperature - tpol)*255/(teq-tpol)
            If red>255
              red=255
            EndIf
            If red<0
              red=0
            EndIf
            blue = 255-red
            green=blue/2
          EndIf
         
        Case #MAPTYPE_SLOPE
          red = geostat(numHmap,i,j)\slope*32
          If red > 255
            red = 255
          EndIf
          green = red:blue = Red ; (greyscale)
         
        Case #MAPTYPE_TERRAIN
          ; Add ice/snow, slope marking
            red=Red(geostat(numHmap,i,j)\colour)
            green=Green(geostat(numHmap,i,j)\colour)
            blue=Blue(geostat(numHmap,i,j)\colour)
         
            ; Ice / snow
            If geostat(numHmap,i,j)\temperature < 0
              If geostat(numHmap,i,j)\temperature < 0
                icegradient= -geostat(numHmap,i,j)\temperature / 5
              EndIf
             
              If icegradient>1
                icegradient=1
              EndIf         
              icecolor=247+Random(8)
             
              red=red*(1-icegradient) + icecolor*icegradient
              green=green*(1-icegradient) + icecolor*icegradient
              blue=blue*(1-icegradient) + icecolor*icegradient       
            EndIf
           
            ; Slopes are darker
            If geostat(numHmap,i,j)\isWater<2   
              red-geostat(numHmap,i,j)\slope * 4
              green-geostat(numHmap,i,j)\slope * 4
              blue-geostat(numHmap,i,j)\slope * 4       
            EndIf     
             
            ; RGB must stay between 0-255
            If red<0
              red=0
            Else
              If red>255
                red=255
              EndIf
            EndIf
            If green<0
              green=0
            Else
              If green>255
                green=255
              EndIf
            EndIf
            If blue<0
              blue=0
            Else
              If blue>255
                blue=255
              EndIf
            EndIf
         
      EndSelect
     
      CompilerIf #USEAPIDRAWING = #True
        color=RGB(red,green,blue)
        *pixel\l=color
        *pixel + 4
      CompilerElse
        color=RGB(blue,green,red)
        Plot(i,j,color)
      CompilerEndIf
    Next i
  Next j
 
  CompilerIf #USEAPIDRAWING = #True
    SetDIBits_(hDC,hBmp,1,gridSize+1,*mem,bmi,#DIB_RGB_COLORS)
  CompilerEndIf
 
  StopDrawing()
 
  CompilerIf #USEAPIDRAWING = #True
    FreeMemory(*mem)
  CompilerEndIf
 
  Debug "Make terrain image: "+Str(ElapsedMilliseconds()-t)+" ms."

  ProcedureReturn numImage
EndProcedure



;- --- Mesh generation ---
;************************************************************************************
; Name: initializeMeshDatas
; Purpose: Creates a list of vertices and a list of triangles for a "plane" mesh
; Parameters:
;   - vertices list
;   - triangles list
;   - sizes x,y,z of the plane
; Return-value: none, but the results are stored in the "vertex" and "triangle" lists()
;************************************************************************************
Procedure.i initializeMeshDatas(List vertex.vertex_struct(),List triangle.triangle_struct(),sizeX.f=1,sizeY.f=1,sizeZ.f=1)

  Protected nbVert.i , nbTri.i              ; Number of vertices and faces
  Protected x.f,y.f,z.f                     ; vertex position
  Protected i.i
 
  ; Read number of vertices and triangles
  Read.i nbVert
  Read.i nbTri
 
  ; Read and store vertices position, normals, uv coords
  For i = 1 To nbVert
    AddElement(vertex())
   
    Read.f x
    Read.f y
    Read.f z
   
    vertex()\x = x * sizex
    vertex()\y = y * sizey
    vertex()\z = z * sizez
   
    Read.f vertex()\nx
    Read.f vertex()\ny
    Read.f vertex()\nz
    Read.f vertex()\u
    Read.f vertex()\v
   
    vertex()\colour = $FFFFFF
   
  Next i   
 
  ;Read and store faces infos
  For i=1 To nbTri
    AddElement(triangle())
   
    Read.i triangle()\v1
    Read.i triangle()\v2
    Read.i triangle()\v3
   
    ; store vertices' pointers, for faster access
    SelectElement(vertex(),triangle()\v1)
    triangle()\ptrV1 = @vertex()
    SelectElement(vertex(),triangle()\v2)
    triangle()\ptrV2 = @vertex()
    SelectElement(vertex(),triangle()\v3)
    triangle()\ptrV3 = @vertex()
  Next i
EndProcedure

;************************************************************************************
; Name: CopyAndRotateVertexList
; Purpose: Copy and rotate a vertex list to make all the sides of the planet from
;          the "top" one. (faster than tesselate the six sides)
; Parameters:
;   - vertices list for "top" side
;   - vertices list to store the result
;   - side you want to create
; Return-value: none, but the result is stored in the "vertex" list()
;************************************************************************************
Procedure copyAndRotateVertexList(List vertex.vertex_struct(),List newVertex.vertex_struct(),side.i) 
  ; Note: This procedure assumes the "vertex" list is the "TOP" one
  CopyList(vertex(),newVertex())
  Select side
    Case #FRONTSIDE
      ForEach newVertex()
        Swap newVertex()\y,newVertex()\z
        newVertex()\y = -newVertex()\y       
      Next newVertex()

    Case #BACKSIDE
      ForEach newVertex()
        Swap newVertex()\y,newVertex()\z
        newVertex()\y = -newVertex()\y
       
        newVertex()\x = -newVertex()\x
        newVertex()\z = -newVertex()\z
      Next newVertex()

    Case #LEFTSIDE
      ForEach newVertex()
        Swap newVertex()\y,newVertex()\z
        newVertex()\y = -newVertex()\y       
       
        Swap newVertex()\x,newVertex()\z
        newVertex()\x = -newVertex()\x
      Next newVertex()

    Case #RIGHTSIDE
      ForEach newVertex()
        Swap newVertex()\y,newVertex()\z
        newVertex()\y = -newVertex()\y       
       
        Swap newVertex()\x,newVertex()\z
        newVertex()\z = -newVertex()\z
      Next newVertex()

    Case #BOTTOMSIDE
      ForEach newVertex()
        newVertex()\y = -newVertex()\y       
        newVertex()\z = -newVertex()\z       
      Next newVertex()
     
  EndSelect
 
EndProcedure

;************************************************************************************
; Name: tesselateMesh
; Purpose: Break each triangle of the mesh into four smaller triangle to increase LOD
; Parameters:
;   - vertices list of the mesh
;   - triangles list of the mesh
;   - number of times you want to repeat the operation
; Return-value: none, but the results are stored in the "vertex" and "triangle" lists()
;************************************************************************************
Procedure.i tesselateMesh(List vertex.vertex_struct(),List triangle.triangle_struct(),depth.i)

  Protected i.i,j.i
  Protected t0.i
 
  Protected NewList newVertex.vertex_struct()
  Protected NewList newTriangle.triangle_struct()
 
  Protected v1_2.vertex_struct,v2_3.vertex_struct,v1_3.vertex_struct
  Protected vertNum1.i,vertNum2.i,vertNum3.i,vertNum1_2.i,vertNum1_3.i,vertNum2_3.i
 
  Protected NewMap antidup.i() ; ugly fix to avoid duplicate vertices...
 
 
  ; Tesselate!
  t0 = ElapsedMilliseconds()
  For i=1 To depth
    ClearList(newVertex())
    ClearList(newTriangle())
    ClearMap(antidup())
   
    ForEach triangle()
      ; Compute new points (in the middle of each side of the triangle)
      MIDDLEPOINT(v1_2,triangle()\ptrV1,triangle()\ptrV2)
      MIDDLEPOINT(v1_3,triangle()\ptrV1,triangle()\ptrV3)
      MIDDLEPOINT(v2_3,triangle()\ptrV2,triangle()\ptrV3)
     
      ; store them in the "new vertex" list, and check if the already exists
      CREATENEWVERTEX(triangle()\ptrV1,vertNum1)
      CREATENEWVERTEX(triangle()\ptrV2,vertNum2)
      CREATENEWVERTEX(triangle()\ptrV3,vertNum3)
      CREATENEWVERTEX(v1_2,vertNum1_2)
      CREATENEWVERTEX(v1_3,vertNum1_3)
      CREATENEWVERTEX(v2_3,vertNum2_3)
     
      ; Create the new triangles
      AddElement(newTriangle())
      newTriangle()\v1 = vertNum1
      newTriangle()\v2 = vertNum1_2
      newTriangle()\v3 = vertNum1_3

      AddElement(newTriangle())
      newTriangle()\v1 = vertNum1_2
      newTriangle()\v2 = vertNum2
      newTriangle()\v3 = vertNum2_3
     
      AddElement(newTriangle())
      newTriangle()\v1 = vertNum1_3
      newTriangle()\v2 = vertNum2_3
      newTriangle()\v3 = vertNum3

      AddElement(newTriangle())
      newTriangle()\v1 = vertNum2_3
      newTriangle()\v2 = vertNum1_3
      newTriangle()\v3 = vertNum1_2
    Next triangle()
   
    ; The new list replaces the old one
    ClearList(vertex())
    ClearList(triangle())
   
    CopyList(newVertex(),vertex())
    CopyList(newTriangle(),triangle())
   
    ; Store vertices' pointers, for faster access
    ForEach triangle()
      SelectElement(vertex(),triangle()\v1)
      triangle()\ptrV1 = @vertex()
      SelectElement(vertex(),triangle()\v2)
      triangle()\ptrV2 = @vertex()
      SelectElement(vertex(),triangle()\v3)
      triangle()\ptrV3 = @vertex()     
    Next triangle()   
   
  Next i
 
  FreeMap(antidup())
 
  Debug "Depth = " + Str(depth)
  Debug "Nb vertices = " + Str(ListSize(vertex()))
  Debug "Nb triangles = " + Str(ListSize(triangle()))
  Debug "Mesh tesselation: "+Str(ElapsedMilliseconds()-t0) + " ms."
 
EndProcedure


;************************************************************************************
; Name: applyHeightmap
; Purpose: Apply a heightmap to a side of the planet
; Parameters:
;   - Heightmap number
;   - vertices list of the mesh
;   - planet size
; Return-value: none, but the result is stored in the "vertex" list()
;************************************************************************************
Procedure applyHeightmap(numHmap.i,List vertex.vertex_struct(),planetRadius.f,reliefExageration.f = 0.2,gridSize.i = #GRIDSIZE,waterLevel.i=#WATERLEVEL)
  Protected height.i
  Protected displacement.f
 
  gridSize-1
  ForEach vertex()
    height = heightmap(numHmap,Int(vertex()\u * gridSize),Int(vertex()\v * gridSize)) - waterLevel
    If height > 0
      displacement = height*reliefExageration
    Else
      displacement = 0
    EndIf
   
    ; This transforms the cube into a sphere
    VEC3_NORMALIZETOLENGTH(vertex(),planetRadius+displacement)
   
  Next vertex()
 
EndProcedure

;************************************************************************************
; Name: weldSideBorders
; Purpose: Adjust the vertices' coords on the borders of the planet "sides", to avoid
;          gaps between the submeshes
; Parameters:
;   none
; Return-value: none, but the results are stored in the "vertex" and "triangle" lists()
;************************************************************************************
Procedure weldSideBorders()
  Protected i.i
 
  ; First , Find and store the vertices which are along the borders of the mesh
  For i=0 To 5
    ForEach planetSide(i)\vertex()     
      If planetSide(i)\vertex()\v = 0
        AddElement(planetSide(i)\topBorder())
        planetSide(i)\topBorder()\ptrVert = @planetSide(i)\vertex()
        planetSide(i)\topBorder()\u = planetSide(i)\vertex()\u
        planetSide(i)\topBorder()\v = planetSide(i)\vertex()\v
      EndIf
      If planetSide(i)\vertex()\v = 1
        AddElement(planetSide(i)\bottomBorder())
        planetSide(i)\bottomBorder()\ptrVert = @planetSide(i)\vertex()
        planetSide(i)\bottomBorder()\u = planetSide(i)\vertex()\u
        planetSide(i)\bottomBorder()\v = planetSide(i)\vertex()\v
      EndIf
      If planetSide(i)\vertex()\u = 0
        AddElement(planetSide(i)\leftBorder())
        planetSide(i)\leftBorder()\ptrVert = @planetSide(i)\vertex()
        planetSide(i)\leftBorder()\u = planetSide(i)\vertex()\u
        planetSide(i)\leftBorder()\v = planetSide(i)\vertex()\v
      EndIf
      If planetSide(i)\vertex()\u = 1
        AddElement(planetSide(i)\rightBorder())
        planetSide(i)\rightBorder()\ptrVert = @planetSide(i)\vertex()
        planetSide(i)\rightBorder()\u = planetSide(i)\vertex()\u
        planetSide(i)\rightBorder()\v = planetSide(i)\vertex()\v
      EndIf
    Next planetSide(i)\vertex()
  Next i
   
  ; Weld top to front/back/left/right
  WELDBORDERVERTICES(#TOPSIDE,bottom,u,#PB_Sort_Ascending,#FRONTSIDE,top,u)
  WELDBORDERVERTICES(#TOPSIDE,left,v,#PB_Sort_Ascending,#LEFTSIDE,top,u)
  WELDBORDERVERTICES(#TOPSIDE,right,v,#PB_Sort_Descending,#RIGHTSIDE,top,u)
  WELDBORDERVERTICES(#TOPSIDE,top,u,#PB_Sort_Descending,#BACKSIDE,top,u)
 
  ; Weld bottom to front/back/left/right
  WELDBORDERVERTICES(#BOTTOMSIDE,top,u,#PB_Sort_Ascending,#FRONTSIDE,bottom,u)
  WELDBORDERVERTICES(#BOTTOMSIDE,left,v,#PB_Sort_Descending,#LEFTSIDE,bottom,u)
  WELDBORDERVERTICES(#BOTTOMSIDE,right,v,#PB_Sort_Ascending,#RIGHTSIDE,bottom,u)
  WELDBORDERVERTICES(#BOTTOMSIDE,bottom,u,#PB_Sort_Descending,#BACKSIDE,bottom,u)
 
  ; Weld front/back/left/right
  WELDBORDERVERTICES(#FRONTSIDE,right,v,#PB_Sort_Ascending,#RIGHTSIDE,left,v)
  WELDBORDERVERTICES(#RIGHTSIDE,right,v,#PB_Sort_Ascending,#BACKSIDE,left,v)
  WELDBORDERVERTICES(#BACKSIDE,right,v,#PB_Sort_Ascending,#LEFTSIDE,left,v)
  WELDBORDERVERTICES(#LEFTSIDE,right,v,#PB_Sort_Ascending,#FRONTSIDE,left,v)
 
  ; Clear the lists to save memory
  For i=0 To 5
    ClearList(planetSide(i)\topBorder())
    ClearList(planetSide(i)\bottomBorder())
    ClearList(planetSide(i)\leftBorder())
    ClearList(planetSide(i)\rightBorder())
  Next i
 
 
EndProcedure


;************************************************************************************
; Name: createMeshFromLists
; Purpose: Create a mesh from the six "vertex()" and "triangle()" lists of the planet.
;          (There's one submesh for each side of the planet)
; Parameters:
;   none
; Return-value: Planet mesh number
;************************************************************************************
Procedure.i createMeshFromLists(Array planetSide.planetSide_struct(1))
  Protected numMesh.i,i.i
 
  ; Create the mesh
  numMesh = CreateMesh(#PB_Any)
 
  ; Add vertices to mesh
  For i =0 To 5
    AddSubMesh()
   
    ForEach planetSide(i)\vertex()   
      MeshVertexPosition(planetSide(i)\vertex()\x,planetSide(i)\vertex()\y,planetSide(i)\vertex()\z)
      ;MeshVertexColor(vertex()\colour)
      MeshVertexNormal(planetSide(i)\vertex()\nx,planetSide(i)\vertex()\ny,planetSide(i)\vertex()\nz)
      MeshVertexTextureCoordinate(planetSide(i)\vertex()\u,planetSide(i)\vertex()\v)
    Next planetSide(i)\vertex()
   
    ; Add triangles to mesh
    ForEach planetSide(i)\triangle()
      MeshFace(planetSide(i)\triangle()\v1,planetSide(i)\triangle()\v2,planetSide(i)\triangle()\v3)
    Next planetSide(i)\triangle()       
   
  Next i
 
  ; Finalize mesh
  FinishMesh(#True)
  NormalizeMesh(numMesh) ; NB: it doesn't normalize between submeshes => visible seams
 
  ; Default material = planet sides textures
  For i =0 To 5
   SetMeshMaterial(numMesh, MaterialID(planetSide(i)\numMaterial),i)
  Next i
 
  ProcedureReturn numMesh
 
EndProcedure   

;************************************************************************************
; Name: createPlanet
; Purpose: Create a planet entity from the "planetSide" array, which contains heighmaps
;          and textures for the six "sides" of the planet  , then attach it to a node   
; Parameters:
;   - "PlanetSide" array
;   - Parent node number
;   - Tesselation level
; Return-value: Planet mesh number
;************************************************************************************
Procedure createPlanet(Array planetSide.planetSide_struct(1), tesselationLevel.i)
  Protected i.i
  Protected planetMesh.i,planetEntity.i
 
 
 ; Clear vertices and triangles list
  For i=0 To 5
    ClearList(planetSide(i)\vertex())
    ClearList(planetSide(i)\triangle())
  Next i

 
  ; Read datas to create the "top" side and tesselate it
  Restore meshDatas
  InitializeMeshDatas(planetSide(#TOPSIDE)\vertex(),planetSide(#TOPSIDE)\triangle(),500,500,500)
  tesselateMesh(planetSide(#TOPSIDE)\vertex(),planetSide(#TOPSIDE)\triangle(),tesselationLevel)
 
  ; Copy it to the other sides
  For i=1 To 5
    CopyAndRotateVertexList(planetSide(#TOPSIDE)\vertex(),planetSide(i)\vertex(),i)
    CopyList(planetSide(#TOPSIDE)\triangle(),planetSide(i)\triangle())
  Next i
 
  ; Apply heightmaps
  For i=0 To 5
    applyHeightmap(i,planetSide(i)\vertex(),500,0.4,#GRIDSIZE,#WATERLEVEL)
  Next i
 
  ; Ensures each side is correctly welded to the next ones
  weldSideBorders()
 
  ; Create the mesh and the entity
  planetMesh = CreateMeshFromLists(planetSide())
  planetEntity = CreateEntity(#PB_Any,MeshID(planetMesh),#PB_Material_None)
  FreeMesh(planetMesh)
 
  ProcedureReturn planetEntity 
EndProcedure


DisableExplicit
;************************************************************************************

;- ---------- Main program ----------
;- Screen initialisation
If InitEngine3D() = 0
   MessageRequester( "Error" , "Can't initialize 3D, check if engine3D.dll is available" , 0 )
   End
EndIf
If InitSprite() = 0 Or InitKeyboard() = 0  Or InitMouse() = 0
   MessageRequester( "Error" , "Can't find DirectX 7.0 or above" , 0 )
   End
EndIf

OpenWindow(0,0, 0, 800 , 600 ,"Telluric Planet Generator v3",#PB_Window_ScreenCentered)
OpenWindowedScreen(WindowID(0),0,0, 800 , 600,0,0,0,#PB_Screen_SmartSynchronization)
KeyboardMode(#PB_Keyboard_International)



;- Heightmaps
; Prepare six heightmaps (one for each side)
dispersion=#GRIDSIZE
t0 = ElapsedMilliseconds()
For i=0 To 5
  ; Compute a "Perlin-like" heightmap
  makeDiamondSquare(#GRIDSIZE,dispersion,#PLANETSEED+i)
 
  ; Voronoi can add a nice "fractured" touch to heightmaps, but can be long to generate
  If #USEVORONOI
    makeVoronoi(#GRIDSIZE,3,3,#PLANETSEED+i)
  EndIf
 
  ; Pass #true as last parameter to use Voronoi
  makeHeightmap(i,#GRIDSIZE,#True,#USEVORONOI) 
Next i
Debug "Heightmaps generation: "+Str(ElapsedMilliseconds()-t0)+" ms."

;- Blend seams between heightmaps
eraseHeightmapSeam(#TOPSIDE,#HEIGHTMAP_BOTTOM,#FRONTSIDE,#HEIGHTMAP_TOP,#GRIDSIZE,25)
eraseHeightmapSeam(#TOPSIDE,#HEIGHTMAP_TOP,#BACKSIDE,#HEIGHTMAP_TOP,#GRIDSIZE,25)
eraseHeightmapSeam(#TOPSIDE,#HEIGHTMAP_LEFT,#LEFTSIDE,#HEIGHTMAP_TOP,#GRIDSIZE,25)
eraseHeightmapSeam(#TOPSIDE,#HEIGHTMAP_RIGHT,#RIGHTSIDE,#HEIGHTMAP_TOP,#GRIDSIZE,25)

eraseHeightmapSeam(#FRONTSIDE,#HEIGHTMAP_LEFT,#LEFTSIDE,#HEIGHTMAP_RIGHT,#GRIDSIZE,25)
eraseHeightmapSeam(#FRONTSIDE,#HEIGHTMAP_RIGHT,#RIGHTSIDE,#HEIGHTMAP_LEFT,#GRIDSIZE,25)
eraseHeightmapSeam(#BACKSIDE,#HEIGHTMAP_RIGHT,#LEFTSIDE,#HEIGHTMAP_LEFT,#GRIDSIZE,25)
eraseHeightmapSeam(#BACKSIDE,#HEIGHTMAP_LEFT,#RIGHTSIDE,#HEIGHTMAP_RIGHT,#GRIDSIZE,25)
 
eraseHeightmapSeam(#BOTTOMSIDE,#HEIGHTMAP_TOP,#FRONTSIDE,#HEIGHTMAP_BOTTOM,#GRIDSIZE,25)
eraseHeightmapSeam(#BOTTOMSIDE,#HEIGHTMAP_BOTTOM,#BACKSIDE,#HEIGHTMAP_BOTTOM,#GRIDSIZE,25)
eraseHeightmapSeam(#BOTTOMSIDE,#HEIGHTMAP_LEFT,#LEFTSIDE,#HEIGHTMAP_BOTTOM,#GRIDSIZE,25)
eraseHeightmapSeam(#BOTTOMSIDE,#HEIGHTMAP_RIGHT,#RIGHTSIDE,#HEIGHTMAP_BOTTOM,#GRIDSIZE,25)


; Read terrain types
i=0
Read.s name
While name <> "XXX END XXX"
  i+1
  terrain(i)\name = name
  Read.f terrain(i)\heightAv
  Read.f terrain(i)\heightVt
  Read.f terrain(i)\slopeAv
  Read.f terrain(i)\slopeVt
  Read.f terrain(i)\temperatureAv
  Read.f terrain(i)\temperatureVt
  Read.i terrain(i)\colour
 
  Read.s name
Wend

;- Textures and materials
Add3DArchive(".",#PB_3DArchive_FileSystem)
For i=0 To 5
   computeGeoStats(i,#GRIDSIZE,#WATERLEVEL,tpol,teq,#PLANETSEED)
   
   planetSide(i)\numTexture = makeTerrainImage(i,#GRIDSIZE,#MAPTYPE_TERRAIN)
   SaveImage(planetSide(i)\numTexture,"temp"+Str(i)+".bmp")
   FreeImage(planetSide(i)\numTexture)
 
  planetSide(i)\numTexture = LoadTexture(#PB_Any,"temp"+Str(i)+".bmp")
  DeleteFile("temp"+Str(i)+".bmp")
 
  planetSide(i)\numMaterial = CreateMaterial(#PB_Any,TextureID(planetSide(i)\numTexture))           
  MaterialFilteringMode(planetSide(i)\numMaterial, #PB_Material_None)
Next i

;- Mesh and entity
planetEntity = createPlanet(planetSide(),tesselation)

; Attach a node, so you can add other elements later (atmosphere, moons, etc..)
CreateNode(0,0,0,0)
AttachNodeObject(0,EntityID(planetEntity))
 

;- Camera
CreateCamera(0, 0, 0 , 100 , 100)
MoveCamera(0,0,1000,1500)
CameraLookAt(0,0,0,0)

;- Light
;WorldShadows(#PB_Shadow_Additive)
AmbientColor($222222)
CreateLight(0,$FFFFFF,400,500,1000)


;- Main loop
Repeat
  While WindowEvent():Wend
  Delay(1)
   
  ;- Mouse
  ExamineMouse()
  anglex+MouseDeltaX()
  angley+MouseDeltaY()
  movez = MouseWheel()

  ;- Keyboard
  ExamineKeyboard()
 
  If KeyboardReleased(#PB_Key_W)
    wireMode = 1-wireMode
    If wireMode = 1
      CameraRenderMode(0,#PB_Camera_Wireframe)
    Else
      CameraRenderMode(0,#PB_Camera_Textured)
    EndIf
  EndIf
 
  If KeyboardReleased(#PB_Key_R) And tesselation > 1
    DetachNodeObject(0,EntityID(planetEntity))
    FreeEntity(planetEntity)
    tesselation-1
    Debug tesselation
    planetEntity = createPlanet(planetSide(),tesselation)
    AttachNodeObject(0,EntityID(planetEntity))
  EndIf
 
  If KeyboardReleased(#PB_Key_T) And tesselation < 10
    DetachNodeObject(0,EntityID(planetEntity))
    FreeEntity(planetEntity)
    tesselation+1
    planetEntity = createPlanet(planetSide(),tesselation)
    AttachNodeObject(0,EntityID(planetEntity))
  EndIf
 
  If KeyboardPushed(#PB_Key_Right)
    RotateNode(0,0,-1,0,#PB_Relative)
    autoRotate = #False
  EndIf
  If KeyboardPushed(#PB_Key_Left)
    RotateNode(0,0,1,0,#PB_Relative)
    autoRotate = #False
  EndIf
 
  If KeyboardPushed(#PB_Key_Up)
    RotateNode(0,1,0,0,#PB_Relative)
    autoRotate = #False
  EndIf
  If KeyboardPushed(#PB_Key_Down)
    RotateNode(0,-1,0,0,#PB_Relative)
    autoRotate = #False
  EndIf
 
 
  MoveCamera(0,0,0,movez * -50)
 
  ; show it all
  RenderWorld()
 
  ; Flip buffers to avoid tearing 
  FlipBuffers()
Until KeyboardPushed(#PB_Key_Escape)
End


;************************************************************************************
;-                                 ---- Datas ----
;************************************************************************************
DataSection
  terrainTypes:
  Data.s "Tundra"
  Data.f 1,4,0,2.5,0,1.5
  Data.i $79A700
 
  Data.s "Grass"
  Data.f 1,3,0,2,18,3.5
  Data.i $5AD00F
 
  Data.s "Forest"
  Data.f 1.1,3,5,3,25,2
  Data.i $3EA300
 
  Data.s "Desert"
  Data.f 1.3,5,0,1.25,50,2
  Data.i $E5DFAD
 
  Data.s "Mountain"
  Data.f 1.7,2,20,3,10,2
  Data.i $BBBBBB
 
  Data.s "XXX END XXX"
 
  meshDatas:
  ; Nb vertices, nb triangles
  Data.i 4,2
  ; Vertices: pos / normals / uv
  Data.f -0.5,0.5,-0.5
  Data.f 0,1,0
  Data.f 0,0
  Data.f 0.5,0.5,-0.5
  Data.f 0,1,0
  Data.f 1,0
  Data.f 0.5,0.5,0.5
  Data.f 0,1,0
  Data.f 1,1
  Data.f -0.5,0.5,0.5
  Data.f 0,1,0
  Data.f 0,1
  ; Faces
  Data.i 2,1,0
  Data.i 0,3,2

EndDataSection

GallyHC · Message par **GallyHC** » ven. 05/juil./2013 11:03

Vraiment sympa GG va avoir de quoi travaille

Cordialement,
GallyHC

SPH · Message par **SPH** » ven. 05/juil./2013 11:04

Assembleur error avec la version 5.20b3

SPH · Message par **SPH** » ven. 05/juil./2013 13:29

Ca marche sur la beta 5
C'est genial

comtois · Message par **comtois** » ven. 05/juil./2013 18:38

Très belle démo.

PureBasic

Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes

Re: Générateur de planètes