src/WeatherCM/FGVoronoi.cpp

   1 /*****************************************************************************
   2
   3  Module:       FGVoronoi.cpp
   4  Author:       Christian Mayer
   5  Date started: 28.05.99
   6  Called by:    main program
   7
   8  ---------- Copyright (C) 1999  Christian Mayer (vader@t-online.de) ----------
   9
  10  This program is free software; you can redistribute it and/or modify it under
  11  the terms of the GNU General Public License as published by the Free Software
  12  Foundation; either version 2 of the License, or (at your option) any later
  13  version.
  14
  15  This program is distributed in the hope that it will be useful, but WITHOUT
  16  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
  17  FOR A PARTICULAR PURPOSE.  See the GNU General Public License for more
  18  details.
  19
  20  You should have received a copy of the GNU General Public License along with
  21  this program; if not, write to the Free Software Foundation, Inc., 59 Temple
  22  Place - Suite 330, Boston, MA  02111-1307, USA.
  23
  24  Further information about the GNU General Public License can also be found on
  25  the world wide web at http://www.gnu.org.
  26
  27 FUNCTIONAL DESCRIPTION
  28 ------------------------------------------------------------------------------
  29 library for Voronoi Diagram calculation based on Steven Fortune 'Sweep2'
  30 FGVoronoi is the wraper to feed the voronoi calulation with a vetor of points
  31 and any class you want, as it uses templates
  32 NOTE: Sweep2 didn't free *any* memory. So I'm doing what I can, but that's not
  33       good enough...
  34
  35 HISTORY
  36 ------------------------------------------------------------------------------
  37 30.05.99   Christian Mayer      Created
  38 16.06.99   Durk Talsma          Portability for Linux
  39 *****************************************************************************/
  40
  41 /****************************************************************************/
  42 /* INCLUDES                                                                 */
  43 /****************************************************************************/
  44 #include "FGVoronoi.h"
  45 #include <Voronoi/voronoi.h>
  46 #include <Voronoi/my_memory.h>
  47
  48 extern "C" {
  49 #include <Voronoi/defs.h>
  50
  51 //forward definitions
  52 void voronoi(int triangulate, struct Site *(*nextsite)());
  53 void geominit(void);
  54 void freeinit(struct Freelist *fl, int size);
  55 struct Site *nextone(void);
  56 bool readsites(PointList input);
  57 };
  58
  59 /****************************************************************************/
  60 /********************************** CODE ************************************/
  61 /****************************************************************************/
  62 FGVoronoiOutputList Voronoiate(const FGVoronoiInputList& input)
  63 {
  64     FGVoronoiOutputList ret_list;
  65
  66     PointList p2ds;
  67
  68     FGVoronoiInputList::iterator it1;
  69
  70     //get the points
  71     for (it1 = (FGVoronoiInputList::iterator)input.begin(); it1 != input.end(); it1++)
  72     {
  73         p2ds.push_back(VoronoiPoint(it1->position.x(), it1->position.y()));
  74     }
  75
  76     cl.clear(); //make sure that it's empty
  77
  78     if (readsites(p2ds) == false)
  79         return ret_list;
  80
  81     freeinit(&sfl, sizeof *sites);
  82     siteidx = 0;
  83     geominit();
  84     voronoi(triangulate, nextone);
  85
  86     _my_free();
  87     /*free(sites);    //to prevent a *big* memory leak...
  88     free(ELhash);   //Sweep2 didn't free *any* memory...
  89     free(PQhash);*/
  90
  91     for (cellList::iterator it2 = cl.begin(); it2 != cl.end(); it2++)
  92     {
  93         it2->sort();
  94         it2->strip();
  95
  96         //uncomment for debugging
  97         //cout << *it2;
  98
  99         Point2DList boundary;
 100
 101         //copy points
 102         PointList::iterator it3 = it2->boundary.begin();
 103
 104         if (it3->infinity == false)
 105             boundary.push_back( *it3 );
 106         else
 107         {
 108             Point2D direction_vector = *it3;
 109             it3++;
 110             boundary.push_back( (*it3) + direction_vector);
 111         }
 112
 113         for (; it3 != it2->boundary.end(); it3++)
 114         {
 115             boundary.push_back( *it3 );
 116
 117         }
 118
 119         it3--;
 120         if (it3->infinity == true)
 121         {
 122             Point2D direction_vector = *it3;
 123             it3--;
 124             Point2D value = *it3;
 125             boundary.pop_back();
 126             boundary.push_back(value + direction_vector);
 127         }
 128
 129         ret_list.push_back(FGVoronoiOutput(boundary, input[it2->ID].value));
 130
 131     }
 132
 133     return ret_list;
 134
 135 }
 136
 137 extern "C"
 138 {
 139
 140 /* return a single in-storage site */
 141 struct Site *nextone(void)
 142 {
 143     struct Site *s;
 144     if(siteidx < nsites)
 145     {
 146         s = &sites[siteidx];
 147         siteidx ++;
 148         return(s);
 149     }
 150     else
 151         return( (struct Site *)NULL);
 152 }
 153
 154
 155 /* sort sites on y, then x, coord */
 156 //int scomp(const struct Point *s1, const struct Point *s2)
 157 int scomp(const void *s1, const void *s2)
 158 {
 159     if(((Point*)s1) -> y < ((Point*)s2) -> y) return(-1);
 160     if(((Point*)s1) -> y > ((Point*)s2) -> y) return(1);
 161     if(((Point*)s1) -> x < ((Point*)s2) -> x) return(-1);
 162     if(((Point*)s1) -> x > ((Point*)s2) -> x) return(1);
 163     return(0);
 164 }
 165
 166
 167
 168 /* read all sites, sort, and compute xmin, xmax, ymin, ymax */
 169 bool readsites(PointList input)
 170 {
 171     int i;
 172
 173     if (input.size() == 0)
 174         return false;   //empty array
 175
 176     PointList::iterator It = input.begin();
 177
 178     nsites=0;
 179     sites = (struct Site *) myalloc(4000*sizeof(*sites));
 180
 181     while(It != input.end())
 182     {
 183         sites[nsites].coord.x = It->x();
 184         sites[nsites].coord.y = It->y();
 185
 186         sites[nsites].sitenbr = nsites;
 187         sites[nsites].refcnt = 0;
 188         nsites ++;
 189         It++;
 190
 191         if (nsites % 4000 == 0)
 192             sites = (struct Site *) my_realloc(sites,(nsites+4000)*sizeof(*sites));
 193     };
 194
 195     qsort(sites, nsites, sizeof (*sites), scomp);
 196     xmin=sites[0].coord.x;
 197     xmax=sites[0].coord.x;
 198
 199     for(i=1; i<nsites; i+=1)
 200     {
 201         if(sites[i].coord.x < xmin)
 202             xmin = sites[i].coord.x;
 203         if(sites[i].coord.x > xmax)
 204             xmax = sites[i].coord.x;
 205     };
 206
 207     ymin = sites[0].coord.y;
 208     ymax = sites[nsites-1].coord.y;
 209
 210     return true;
 211 }
 212
 213
 214 }
 215
 216