PLplot 5.9.6
|
00001 /* $Id$ 00002 * 00003 * Continental Outline and Political Boundary Backgrounds 00004 * 00005 * Some plots need a geographical background such as the global 00006 * surface temperatures or the population density. The routine 00007 * plmap() will draw one of the following backgrounds: continental 00008 * outlines, political boundaries, the United States, and the United 00009 * States with the continental outlines. The routine plmeridians() 00010 * will add the latitudes and longitudes to the background. After 00011 * the background has been drawn, one can use a contour routine or a 00012 * symbol plotter to finish off the plot. 00013 * 00014 * Copyright (C) 1991, 1993, 1994 Wesley Ebisuzaki 00015 * Copyright (C) 1994, 2000, 2001 Maurice LeBrun 00016 * Copyright (C) 1999 Geoffrey Furnish 00017 * Copyright (C) 2000, 2001, 2002 Alan W. Irwin 00018 * Copyright (C) 2001 Andrew Roach 00019 * Copyright (C) 2001, 2004 Rafael Laboissiere 00020 * Copyright (C) 2002 Vincent Darley 00021 * Copyright (C) 2003 Joao Cardoso 00022 * 00023 * This file is part of PLplot. 00024 * 00025 * PLplot is free software; you can redistribute it and/or modify 00026 * it under the terms of the GNU Library General Public License 00027 * as published by the Free Software Foundation; version 2 of the 00028 * License. 00029 * 00030 * PLplot is distributed in the hope that it will be useful, but 00031 * WITHOUT ANY WARRANTY; without even the implied warranty of 00032 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00033 * GNU Library General Public License for more details. 00034 * 00035 * You should have received a copy of the GNU Library General Public 00036 * License along with this library; if not, write to the Free Software 00037 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 00038 * USA 00039 * 00040 */ 00041 00042 #include "plplotP.h" 00043 00044 /*----------------------------------------------------------------------*\ 00045 * void plmap(void (*mapform)(PLINT, PLFLT *, PLFLT *), const char *type, 00046 * PLFLT minlong, PLFLT maxlong, PLFLT minlat, PLFLT maxlat); 00047 * 00048 * plot continental outline in world coordinates 00049 * 00050 * v1.4: machine independant version 00051 * v1.3: replaced plcontinent by plmap, added plmeridians 00052 * v1.2: 2 arguments: mapform, type of plot 00053 * 00054 * mapform(PLINT n, PLFLT *x, PLFLT *y) is a routine to transform the 00055 * coordinate longitudes and latitudes to a plot coordinate system. By 00056 * using this transform, we can change from a longitude, latitude 00057 * coordinate to a polar stereographic project, for example. Initially, 00058 * x[0]..[n-1] are the longitudes and y[0]..y[n-1] are the corresponding 00059 * latitudes. After the call to mapform(), x[] and y[] should be replaced 00060 * by the corresponding plot coordinates. If no transform is desired, 00061 * mapform can be replaced by NULL. 00062 * 00063 * type is a character string. The value of this parameter determines the 00064 * type of background. The possible values are, 00065 * 00066 * "globe" continental outlines 00067 * "usa" USA and state boundaries 00068 * "cglobe" continental outlines and countries 00069 * "usaglobe" USA, state boundaries and continental outlines 00070 * 00071 * minlong, maxlong are the values of the longitude on the left and right 00072 * side of the plot, respectively. The value of minlong must be less than 00073 * the values of maxlong, and the values of maxlong-minlong must be less 00074 * or equal to 360. 00075 * 00076 * minlat, maxlat are the minimum and maximum latitudes to be plotted on 00077 * the background. One can always use -90.0 and 90.0 as the boundary 00078 * outside the plot window will be automatically eliminated. However, the 00079 * program will be faster if one can reduce the size of the background 00080 * plotted. 00081 \*----------------------------------------------------------------------*/ 00082 00083 #define MAP_FILE ".map" 00084 #define OFFSET ( 180 * 100 ) 00085 #define SCALE 100.0 00086 #define W_BUFSIZ ( 32 * 1024 ) 00087 00088 void 00089 plmap( void ( *mapform )( PLINT, PLFLT *, PLFLT * ), const char *type, 00090 PLFLT minlong, PLFLT maxlong, PLFLT minlat, PLFLT maxlat ) 00091 { 00092 PLINT wrap, sign; 00093 int i, j; 00094 PLFLT bufx[200], bufy[200], x[2], y[2]; 00095 short int test[400]; 00096 register PDFstrm *in; 00097 char filename[100]; 00098 00099 unsigned char n_buff[2], buff[800]; 00100 int n; 00101 long int t; 00102 00103 /* 00104 * read map outline 00105 */ 00106 strncpy( filename, type, 99 ); 00107 filename[99] = '\0'; 00108 strncat( filename, MAP_FILE, 99 - strlen( filename ) ); 00109 00110 if ( ( in = plLibOpenPdfstrm( filename ) ) == NULL ) 00111 return; 00112 00113 for (;; ) 00114 { 00115 /* read in # points in segment */ 00116 if ( pdf_rdx( n_buff, sizeof ( unsigned char ) * 2, in ) == 0 ) 00117 break; 00118 n = ( n_buff[0] << 8 ) + n_buff[1]; 00119 if ( n == 0 ) 00120 break; 00121 00122 pdf_rdx( buff, sizeof ( unsigned char ) * 4 * n, in ); 00123 if ( n == 1 ) 00124 continue; 00125 00126 for ( j = i = 0; i < n; i++, j += 2 ) 00127 { 00128 t = ( buff[j] << 8 ) + buff[j + 1]; 00129 bufx[i] = ( t - OFFSET ) / SCALE; 00130 } 00131 for ( i = 0; i < n; i++, j += 2 ) 00132 { 00133 t = ( buff[j] << 8 ) + buff[j + 1]; 00134 bufy[i] = ( t - OFFSET ) / SCALE; 00135 } 00136 00137 for ( i = 0; i < n; i++ ) 00138 { 00139 while ( bufx[i] < minlong ) 00140 { 00141 bufx[i] += 360.0; 00142 } 00143 while ( bufx[i] > maxlong ) 00144 { 00145 bufx[i] -= 360.0; 00146 } 00147 } 00148 00149 /* remove last 2 points if both outside of domain and won't plot */ 00150 00151 /* AR: 18/11/01 00152 * I have commented out the next section which supposedly 00153 * removes points that do not plot within the domain. 00154 * 00155 * This code appears at any rate to be superseded by the next 00156 * block of code that checks for wrapping problems. Removing 00157 * this code seems to have fixed up the problems with mapping 00158 * function, but I do not wish to delete it outright just now in 00159 * case I have managed to miss something. 00160 */ 00161 00162 /* while (n > 1) { 00163 * if ((bufx[n-1] < minlong && bufx[n-2] < minlong) || 00164 * (bufx[n-1] > maxlong && bufx[n-2] > maxlong) || 00165 * (bufy[n-1] < minlat && bufy[n-2] < minlat) || 00166 * (bufy[n-1] > maxlat && bufy[n-2] > maxlat)) { 00167 * n--; 00168 * } 00169 * else { 00170 * break; 00171 * } 00172 * } 00173 * if (n <= 1) continue; 00174 */ 00175 00176 if ( mapform != NULL ) 00177 ( *mapform )( n, bufx, bufy ); /* moved transformation to here */ 00178 /* so bound-checking worked right */ 00179 00180 wrap = 0; 00181 /* check for wrap around problems */ 00182 for ( i = 0; i < n - 1; i++ ) 00183 { 00184 /* jc: this code is wrong, as the bufx/y are floats that are 00185 * converted to ints before abs() is called. Thus, small 00186 * differences are masked off. The proof that it is wrong is 00187 * that when replacing abs() with fabs(), as it should be, 00188 * the code works the wrong way. What a proof :-) 00189 * 00190 * test[i] = abs((bufx[i]-bufx[i+1])) > abs(bufy[i]/3); 00191 * 00192 * If the intended behaviour is *really* that, than an 00193 * explicit cast should be used to help other programmers, as 00194 * is done bellow!!! 00195 */ 00196 00197 test[i] = abs( (int) ( bufx[i] - bufx[i + 1] ) ) > abs( (int) bufy[i] / 3 ); /* Changed this from 30 degrees so it is now "polar sensitive" */ 00198 if ( test[i] ) 00199 wrap = 1; 00200 } 00201 00202 if ( wrap == 0 ) 00203 { 00204 plline( n, bufx, bufy ); 00205 } 00206 else 00207 { 00208 for ( i = 0; i < n - 1; i++ ) 00209 { 00210 x[0] = bufx[i]; 00211 x[1] = bufx[i + 1]; 00212 y[0] = bufy[i]; 00213 y[1] = bufy[i + 1]; 00214 if ( test[i] == 0 ) 00215 { 00216 plline( 2, x, y ); 00217 } 00218 else /* this code seems to supercede the block commented out above */ 00219 { 00220 /* segment goes off the edge */ 00221 sign = ( x[1] > x[0] ) ? 1 : -1; 00222 x[1] -= sign * 360.0; 00223 plline( 2, x, y ); 00224 x[0] = bufx[i]; 00225 x[1] = bufx[i + 1]; 00226 y[0] = bufy[i]; 00227 y[1] = bufy[i + 1]; 00228 x[0] += sign * 360.0; 00229 plline( 2, x, y ); 00230 } 00231 } 00232 } 00233 } 00234 /* Close map file */ 00235 pdf_close( in ); 00236 } 00237 00238 /*----------------------------------------------------------------------*\ 00239 * void plmeridians(void (*mapform)(PLINT, PLFLT *, PLFLT *), 00240 * PLFLT dlong, PLFLT dlat, PLFLT minlong, PLFLT maxlong, 00241 * PLFLT minlat, PLFLT maxlat); 00242 * 00243 * Plot the latitudes and longitudes on the background. The lines 00244 * are plotted in the current color and line style. 00245 * 00246 * mapform(PLINT n, PLFLT *x, PLFLT *y) is a routine to transform the 00247 * coordinate longitudes and latitudes to a plot coordinate system. By 00248 * using this transform, we can change from a longitude, latitude 00249 * coordinate to a polar stereographic project, for example. Initially, 00250 * x[0]..x[n-1] are the longitudes and y[0]..y[n-1] are the corresponding 00251 * latitudes. After the call to mapform(), x[] and y[] should be replaced 00252 * by the corresponding plot coordinates. If no transform is desired, 00253 * mapform can be replaced by NULL. 00254 * 00255 * dlat, dlong are the interval in degrees that the latitude and longitude 00256 * lines are to be plotted. 00257 * 00258 * minlong, maxlong are the values of the longitude on the left and right 00259 * side of the plot, respectively. The value of minlong must be less than 00260 * the values of maxlong, and the values of maxlong-minlong must be less 00261 * or equal to 360. 00262 * 00263 * minlat, maxlat are the minimum and maximum latitudes to be plotted on 00264 * the background. One can always use -90.0 and 90.0 as the boundary 00265 * outside the plot window will be automatically eliminated. However, the 00266 * program will be faster if one can reduce the size of the background 00267 * plotted. 00268 \*----------------------------------------------------------------------*/ 00269 00270 #define NSEG 100 00271 00272 void 00273 plmeridians( void ( *mapform )( PLINT, PLFLT *, PLFLT * ), 00274 PLFLT dlong, PLFLT dlat, 00275 PLFLT minlong, PLFLT maxlong, PLFLT minlat, PLFLT maxlat ) 00276 { 00277 PLFLT yy, xx, temp, x[2], y[2], dx, dy; 00278 00279 if ( minlong > maxlong ) 00280 { 00281 temp = minlong; 00282 minlong = maxlong; 00283 maxlong = temp; 00284 } 00285 if ( minlat > maxlat ) 00286 { 00287 temp = minlat; 00288 minlat = maxlat; 00289 maxlat = temp; 00290 } 00291 dx = ( maxlong - minlong ) / NSEG; 00292 dy = ( maxlat - minlat ) / NSEG; 00293 00294 /* latitudes */ 00295 00296 for ( yy = dlat * ceil( minlat / dlat ); yy <= maxlat; yy += dlat ) 00297 { 00298 if ( mapform == NULL ) 00299 { 00300 plpath( NSEG, minlong, yy, maxlong, yy ); 00301 } 00302 else 00303 { 00304 for ( xx = minlong; xx < maxlong; xx += dx ) 00305 { 00306 y[0] = y[1] = yy; 00307 x[0] = xx; 00308 x[1] = xx + dx; 00309 ( *mapform )( 2, x, y ); 00310 plline( 2, x, y ); 00311 } 00312 } 00313 } 00314 00315 /* longitudes */ 00316 00317 for ( xx = dlong * ceil( minlong / dlong ); xx <= maxlong; xx += dlong ) 00318 { 00319 if ( mapform == NULL ) 00320 { 00321 plpath( NSEG, xx, minlat, xx, maxlat ); 00322 } 00323 else 00324 { 00325 for ( yy = minlat; yy < maxlat; yy += dy ) 00326 { 00327 x[0] = x[1] = xx; 00328 y[0] = yy; 00329 y[1] = yy + dy; 00330 ( *mapform )( 2, x, y ); 00331 plline( 2, x, y ); 00332 } 00333 } 00334 } 00335 }