00001 /****************************************************************************** 00002 * $Id: geotiff_proj4_c-source.html 385 2001-03-05 04:58:33Z warmerda $ 00003 * 00004 * Project: libgeotiff 00005 * Purpose: Code to convert a normalized GeoTIFF definition into a PROJ.4 00006 * (OGDI) compatible projection string. 00007 * Author: Frank Warmerdam, warmerda@home.com 00008 * 00009 ****************************************************************************** 00010 * Copyright (c) 1999, Frank Warmerdam 00011 * 00012 * Permission is hereby granted, free of charge, to any person obtaining a 00013 * copy of this software and associated documentation files (the "Software"), 00014 * to deal in the Software without restriction, including without limitation 00015 * the rights to use, copy, modify, merge, publish, distribute, sublicense, 00016 * and/or sell copies of the Software, and to permit persons to whom the 00017 * Software is furnished to do so, subject to the following conditions: 00018 * 00019 * The above copyright notice and this permission notice shall be included 00020 * in all copies or substantial portions of the Software. 00021 * 00022 * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS 00023 * OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, 00024 * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL 00025 * THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER 00026 * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING 00027 * FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER 00028 * DEALINGS IN THE SOFTWARE. 00029 ****************************************************************************** 00030 * 00031 * $Log$ 00031 * Revision 1.2 2001/03/05 04:58:33 warmerda 00031 * updated 00031 * 00032 * Revision 1.16 2000/12/05 19:21:45 warmerda 00033 * added cassini support 00034 * 00035 * Revision 1.15 2000/12/05 17:44:41 warmerda 00036 * Use +R_A for Miller and VanDerGrinten 00037 * 00038 * Revision 1.14 2000/10/13 18:06:51 warmerda 00039 * added econic support for PROJ.4 translation 00040 * 00041 * Revision 1.13 2000/09/15 19:30:48 warmerda 00042 * *** empty log message *** 00043 * 00044 * Revision 1.12 2000/09/15 18:21:07 warmerda 00045 * Fixed order of parameters for LCC 2SP. When parameters 00046 * were read from EPSG CSV files the standard parallels and origin 00047 * were mixed up. This affects alot of state plane zones! 00048 * 00049 * Revision 1.11 2000/06/06 17:39:45 warmerda 00050 * Modify to work with projUV version of library. 00051 * 00052 * Revision 1.10 1999/07/06 15:05:51 warmerda 00053 * Fixed up LCC_1SP notes. 00054 * 00055 * Revision 1.9 1999/05/04 16:24:49 warmerda 00056 * Fixed projection string formating with zones. 00057 * 00058 * Revision 1.8 1999/05/04 12:27:01 geotiff 00059 * only emit proj unsupported warning if DEBUG defined 00060 * 00061 * Revision 1.7 1999/05/04 03:14:59 warmerda 00062 * fixed use of foot instead of ft for units 00063 * 00064 * Revision 1.6 1999/05/03 17:50:31 warmerda 00065 * avoid warnings on IRIX 00066 * 00067 * Revision 1.5 1999/04/29 23:02:24 warmerda 00068 * added mapsys utm test. 00069 * 00070 * Revision 1.4 1999/03/18 21:35:42 geotiff 00071 * Added reprojection functions 00072 * 00073 * Revision 1.3 1999/03/10 18:11:17 geotiff 00074 * Removed comment about this not being the master ... now it is. 00075 * 00076 * Revision 1.2 1999/03/10 18:10:27 geotiff 00077 * Avoid use of cpl_serv.h and CPLStrdup(). 00078 * 00079 * Revision 1.1 1999/03/10 15:20:43 geotiff 00080 * New 00081 * 00082 */ 00083 00084 #include "cpl_serv.h" 00085 #include "geotiff.h" 00086 #include "geo_normalize.h" 00087 #include "geovalues.h" 00088 00089 /************************************************************************/ 00090 /* GTIFGetProj4Defn() */ 00091 /************************************************************************/ 00092 00093 char * GTIFGetProj4Defn( GTIFDefn * psDefn ) 00094 00095 { 00096 char szProjection[512]; 00097 char szUnits[24]; 00098 double dfFalseEasting, dfFalseNorthing; 00099 00100 szProjection[0] = '\0'; 00101 00102 /* ==================================================================== */ 00103 /* Translate the units of measure. */ 00104 /* */ 00105 /* Note that even with a +units, or +to_meter in effect, it is */ 00106 /* still assumed that all the projection parameters are in */ 00107 /* meters. */ 00108 /* ==================================================================== */ 00109 if( psDefn->UOMLength == Linear_Meter ) 00110 { 00111 strcpy( szUnits, "+units=m " ); 00112 } 00113 else if( psDefn->UOMLength == Linear_Foot ) 00114 { 00115 strcpy( szUnits, "+units=ft " ); 00116 } 00117 else if( psDefn->UOMLength == Linear_Foot_US_Survey ) 00118 { 00119 strcpy( szUnits, "+units=us-ft " ); 00120 } 00121 else if( psDefn->UOMLength == Linear_Foot_Indian ) 00122 { 00123 strcpy( szUnits, "+units=ind-ft " ); 00124 } 00125 else if( psDefn->UOMLength == Linear_Link ) 00126 { 00127 strcpy( szUnits, "+units=link " ); 00128 } 00129 else if( psDefn->UOMLength == Linear_Yard_Indian) 00130 { 00131 strcpy( szUnits, "+units=ind-yd " ); 00132 } 00133 else if( psDefn->UOMLength == Linear_Fathom ) 00134 { 00135 strcpy( szUnits, "+units=fath " ); 00136 } 00137 else if( psDefn->UOMLength == Linear_Mile_International_Nautical ) 00138 { 00139 strcpy( szUnits, "+units=kmi " ); 00140 } 00141 else 00142 { 00143 sprintf( szUnits, "+to_meter=%.10f", psDefn->UOMLengthInMeters ); 00144 } 00145 00146 /* -------------------------------------------------------------------- */ 00147 /* false easting and northing are in meters and that is what */ 00148 /* PROJ.4 wants regardless of the linear units. */ 00149 /* -------------------------------------------------------------------- */ 00150 dfFalseEasting = psDefn->ProjParm[5]; 00151 dfFalseNorthing = psDefn->ProjParm[6]; 00152 00153 /* ==================================================================== */ 00154 /* Handle general projection methods. */ 00155 /* ==================================================================== */ 00156 00157 /* -------------------------------------------------------------------- */ 00158 /* UTM - special case override on transverse mercator so things */ 00159 /* will be more meaningful to the user. */ 00160 /* -------------------------------------------------------------------- */ 00161 if( psDefn->MapSys == MapSys_UTM_North ) 00162 { 00163 sprintf( szProjection+strlen(szProjection), 00164 "+proj=utm +zone=%d ", 00165 psDefn->Zone ); 00166 } 00167 00168 /* -------------------------------------------------------------------- */ 00169 /* Transverse Mercator */ 00170 /* -------------------------------------------------------------------- */ 00171 else if( psDefn->CTProjection == CT_TransverseMercator ) 00172 { 00173 sprintf( szProjection+strlen(szProjection), 00174 "+proj=tmerc +lat_0=%.9f +lon_0=%.9f +k=%f +x_0=%.3f +y_0=%.3f ", 00175 psDefn->ProjParm[0], 00176 psDefn->ProjParm[1], 00177 psDefn->ProjParm[4], 00178 dfFalseEasting, 00179 dfFalseNorthing ); 00180 } 00181 00182 /* -------------------------------------------------------------------- */ 00183 /* Mercator */ 00184 /* -------------------------------------------------------------------- */ 00185 else if( psDefn->CTProjection == CT_Mercator ) 00186 { 00187 sprintf( szProjection+strlen(szProjection), 00188 "+proj=merc +lat_ts=%.9f +lon_0=%.9f +k=%f +x_0=%.3f +y_0=%.3f ", 00189 psDefn->ProjParm[0], 00190 psDefn->ProjParm[1], 00191 psDefn->ProjParm[4], 00192 dfFalseEasting, 00193 dfFalseNorthing ); 00194 } 00195 00196 /* -------------------------------------------------------------------- */ 00197 /* Cassini/Soldner */ 00198 /* -------------------------------------------------------------------- */ 00199 else if( psDefn->CTProjection == CT_CassiniSoldner ) 00200 { 00201 sprintf( szProjection+strlen(szProjection), 00202 "+proj=cass +lat_0=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00203 psDefn->ProjParm[0], 00204 psDefn->ProjParm[1], 00205 dfFalseEasting, 00206 dfFalseNorthing ); 00207 } 00208 00209 /* -------------------------------------------------------------------- */ 00210 /* Oblique Stereographic - Should this really map onto */ 00211 /* Stereographic? */ 00212 /* -------------------------------------------------------------------- */ 00213 else if( psDefn->CTProjection == CT_ObliqueStereographic ) 00214 { 00215 sprintf( szProjection+strlen(szProjection), 00216 "+proj=stere +lat_0=%.9f +lon_0=%.9f +k=%f +x_0=%.3f +y_0=%.3f ", 00217 psDefn->ProjParm[0], 00218 psDefn->ProjParm[1], 00219 psDefn->ProjParm[4], 00220 dfFalseEasting, 00221 dfFalseNorthing ); 00222 } 00223 00224 /* -------------------------------------------------------------------- */ 00225 /* Stereographic */ 00226 /* -------------------------------------------------------------------- */ 00227 else if( psDefn->CTProjection == CT_Stereographic ) 00228 { 00229 sprintf( szProjection+strlen(szProjection), 00230 "+proj=stere +lat_0=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00231 psDefn->ProjParm[0], 00232 psDefn->ProjParm[1], 00233 dfFalseEasting, 00234 dfFalseNorthing ); 00235 } 00236 00237 /* -------------------------------------------------------------------- */ 00238 /* Polar Stereographic */ 00239 /* -------------------------------------------------------------------- */ 00240 else if( psDefn->CTProjection == CT_PolarStereographic ) 00241 { 00242 sprintf( szProjection+strlen(szProjection), 00243 "+proj=stere +lat_0=%.9f +lon_0=%.9f +k=%.9f +x_0=%.3f +y_0=%.3f ", 00244 psDefn->ProjParm[0], 00245 psDefn->ProjParm[1], 00246 psDefn->ProjParm[4], 00247 dfFalseEasting, 00248 dfFalseNorthing ); 00249 } 00250 00251 /* -------------------------------------------------------------------- */ 00252 /* Equirectangular */ 00253 /* -------------------------------------------------------------------- */ 00254 else if( psDefn->CTProjection == CT_Equirectangular ) 00255 { 00256 sprintf( szProjection+strlen(szProjection), 00257 "+proj=eqc +lat_ts=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00258 psDefn->ProjParm[0], 00259 psDefn->ProjParm[1], 00260 dfFalseEasting, 00261 dfFalseNorthing ); 00262 } 00263 00264 /* -------------------------------------------------------------------- */ 00265 /* Gnomonic */ 00266 /* -------------------------------------------------------------------- */ 00267 else if( psDefn->CTProjection == CT_Gnomonic ) 00268 { 00269 sprintf( szProjection+strlen(szProjection), 00270 "+proj=gnom +lat_0=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00271 psDefn->ProjParm[0], 00272 psDefn->ProjParm[1], 00273 dfFalseEasting, 00274 dfFalseNorthing ); 00275 } 00276 00277 /* -------------------------------------------------------------------- */ 00278 /* Orthographic */ 00279 /* -------------------------------------------------------------------- */ 00280 else if( psDefn->CTProjection == CT_Orthographic ) 00281 { 00282 sprintf( szProjection+strlen(szProjection), 00283 "+proj=ortho +lat_0=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00284 psDefn->ProjParm[0], 00285 psDefn->ProjParm[1], 00286 dfFalseEasting, 00287 dfFalseNorthing ); 00288 } 00289 00290 /* -------------------------------------------------------------------- */ 00291 /* Lambert Azimuthal Equal Area */ 00292 /* -------------------------------------------------------------------- */ 00293 else if( psDefn->CTProjection == CT_LambertAzimEqualArea ) 00294 { 00295 sprintf( szProjection+strlen(szProjection), 00296 "+proj=laea +lat_0=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00297 psDefn->ProjParm[0], 00298 psDefn->ProjParm[1], 00299 dfFalseEasting, 00300 dfFalseNorthing ); 00301 } 00302 00303 /* -------------------------------------------------------------------- */ 00304 /* Azimuthal Equidistant */ 00305 /* -------------------------------------------------------------------- */ 00306 else if( psDefn->CTProjection == CT_AzimuthalEquidistant ) 00307 { 00308 sprintf( szProjection+strlen(szProjection), 00309 "+proj=aeqd +lat_0=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00310 psDefn->ProjParm[0], 00311 psDefn->ProjParm[1], 00312 dfFalseEasting, 00313 dfFalseNorthing ); 00314 } 00315 00316 /* -------------------------------------------------------------------- */ 00317 /* Miller Cylindrical */ 00318 /* -------------------------------------------------------------------- */ 00319 else if( psDefn->CTProjection == CT_MillerCylindrical ) 00320 { 00321 sprintf( szProjection+strlen(szProjection), 00322 "+proj=mill +lat_0=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f +R_A ", 00323 psDefn->ProjParm[0], 00324 psDefn->ProjParm[1], 00325 dfFalseEasting, 00326 dfFalseNorthing ); 00327 } 00328 00329 /* -------------------------------------------------------------------- */ 00330 /* Polyconic */ 00331 /* -------------------------------------------------------------------- */ 00332 else if( psDefn->CTProjection == CT_Polyconic ) 00333 { 00334 sprintf( szProjection+strlen(szProjection), 00335 "+proj=poly +lat_0=%.9f +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00336 psDefn->ProjParm[0], 00337 psDefn->ProjParm[1], 00338 dfFalseEasting, 00339 dfFalseNorthing ); 00340 } 00341 00342 /* -------------------------------------------------------------------- */ 00343 /* AlbersEqualArea */ 00344 /* -------------------------------------------------------------------- */ 00345 else if( psDefn->CTProjection == CT_AlbersEqualArea ) 00346 { 00347 sprintf( szProjection+strlen(szProjection), 00348 "+proj=aea +lat_1=%.9f +lat_2=%.9f +lat_0=%.9f +lon_0=%.9f" 00349 " +x_0=%.3f +y_0=%.3f ", 00350 psDefn->ProjParm[0], 00351 psDefn->ProjParm[1], 00352 psDefn->ProjParm[2], 00353 psDefn->ProjParm[3], 00354 dfFalseEasting, 00355 dfFalseNorthing ); 00356 } 00357 00358 /* -------------------------------------------------------------------- */ 00359 /* EquidistantConic */ 00360 /* -------------------------------------------------------------------- */ 00361 else if( psDefn->CTProjection == CT_EquidistantConic ) 00362 { 00363 sprintf( szProjection+strlen(szProjection), 00364 "+proj=eqdc +lat_1=%.9f +lat_2=%.9f +lat_0=%.9f +lon_0=%.9f" 00365 " +x_0=%.3f +y_0=%.3f ", 00366 psDefn->ProjParm[0], 00367 psDefn->ProjParm[1], 00368 psDefn->ProjParm[2], 00369 psDefn->ProjParm[3], 00370 dfFalseEasting, 00371 dfFalseNorthing ); 00372 } 00373 00374 /* -------------------------------------------------------------------- */ 00375 /* Robinson */ 00376 /* -------------------------------------------------------------------- */ 00377 else if( psDefn->CTProjection == CT_Robinson ) 00378 { 00379 sprintf( szProjection+strlen(szProjection), 00380 "+proj=robin +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00381 psDefn->ProjParm[1], 00382 dfFalseEasting, 00383 dfFalseNorthing ); 00384 } 00385 00386 /* -------------------------------------------------------------------- */ 00387 /* VanDerGrinten */ 00388 /* -------------------------------------------------------------------- */ 00389 else if( psDefn->CTProjection == CT_VanDerGrinten ) 00390 { 00391 sprintf( szProjection+strlen(szProjection), 00392 "+proj=vandg +lon_0=%.9f +x_0=%.3f +y_0=%.3f +R_A ", 00393 psDefn->ProjParm[1], 00394 dfFalseEasting, 00395 dfFalseNorthing ); 00396 } 00397 00398 /* -------------------------------------------------------------------- */ 00399 /* Sinusoidal */ 00400 /* -------------------------------------------------------------------- */ 00401 else if( psDefn->CTProjection == CT_Sinusoidal ) 00402 { 00403 sprintf( szProjection+strlen(szProjection), 00404 "+proj=sinu +lon_0=%.9f +x_0=%.3f +y_0=%.3f ", 00405 psDefn->ProjParm[1], 00406 dfFalseEasting, 00407 dfFalseNorthing ); 00408 } 00409 00410 /* -------------------------------------------------------------------- */ 00411 /* LambertConfConic_2SP */ 00412 /* -------------------------------------------------------------------- */ 00413 else if( psDefn->CTProjection == CT_LambertConfConic_2SP ) 00414 { 00415 sprintf( szProjection+strlen(szProjection), 00416 "+proj=lcc +lat_0=%.9f +lon_0=%.9f +lat_1=%.9f +lat_2=%.9f " 00417 " +x_0=%.3f +y_0=%.3f ", 00418 psDefn->ProjParm[0], 00419 psDefn->ProjParm[1], 00420 psDefn->ProjParm[2], 00421 psDefn->ProjParm[3], 00422 dfFalseEasting, 00423 dfFalseNorthing ); 00424 } 00425 00426 /* -------------------------------------------------------------------- */ 00427 /* LambertConfConic_1SP */ 00428 /* -------------------------------------------------------------------- */ 00429 else if( psDefn->CTProjection == CT_LambertConfConic_1SP ) 00430 { 00431 /* this appears to be an unsupported formulation with PROJ.4 */ 00432 00433 /* to at least some degree this can be treated similarly to 00434 * the 2SP case. 00435 * 00436 * See http://www.mentorsoftwareinc.com/CC/asknorm/ASK0699.HTM#Q2 00437 */ 00438 } 00439 00440 /* -------------------------------------------------------------------- */ 00441 /* NewZealandMapGrid */ 00442 /* -------------------------------------------------------------------- */ 00443 else if( psDefn->CTProjection == CT_NewZealandMapGrid ) 00444 { 00445 /* this appears to be an unsupported formulation with PROJ.4 */ 00446 } 00447 00448 /* -------------------------------------------------------------------- */ 00449 /* Transverse Mercator - south oriented. */ 00450 /* -------------------------------------------------------------------- */ 00451 else if( psDefn->CTProjection == CT_TransvMercator_SouthOriented ) 00452 { 00453 /* this appears to be an unsupported formulation with PROJ.4 */ 00454 } 00455 00456 /* -------------------------------------------------------------------- */ 00457 /* ObliqueMercator (Hotine) */ 00458 /* -------------------------------------------------------------------- */ 00459 else if( psDefn->CTProjection == CT_ObliqueMercator ) 00460 { 00461 /* not clear how ProjParm[3] - angle from rectified to skewed grid - 00462 should be applied ... see the +not_rot flag for PROJ.4. 00463 Just ignoring for now. */ 00464 00465 sprintf( szProjection+strlen(szProjection), 00466 "+proj=omerc +lat_0=%.9f +lonc=%.9f +alpha=%.9f" 00467 " +k=%.9f +x_0=%.3f +y_0=%.3f ", 00468 psDefn->ProjParm[0], 00469 psDefn->ProjParm[1], 00470 psDefn->ProjParm[2], 00471 psDefn->ProjParm[4], 00472 psDefn->ProjParm[5], 00473 psDefn->ProjParm[6] ); 00474 } 00475 00476 /* ==================================================================== */ 00477 /* Handle ellipsoid information. */ 00478 /* ==================================================================== */ 00479 if( psDefn->Ellipsoid == Ellipse_WGS_84 ) 00480 strcat( szProjection, "+ellps=WGS84 " ); 00481 else if( psDefn->Ellipsoid == Ellipse_Clarke_1866 ) 00482 strcat( szProjection, "+ellps=clrk66 " ); 00483 else if( psDefn->Ellipsoid == Ellipse_Clarke_1880 ) 00484 strcat( szProjection, "+ellps=clrk80 " ); 00485 else if( psDefn->Ellipsoid == Ellipse_GRS_1980 ) 00486 strcat( szProjection, "+ellps=GRS80 " ); 00487 else 00488 { 00489 if( psDefn->SemiMajor != 0.0 && psDefn->SemiMinor != 0.0 ) 00490 { 00491 sprintf( szProjection+strlen(szProjection), 00492 "+a=%.3f +b=%.3f ", 00493 psDefn->SemiMajor, 00494 psDefn->SemiMinor ); 00495 } 00496 } 00497 00498 strcat( szProjection, szUnits ); 00499 00500 return( strdup( szProjection ) ); 00501 } 00502 00503 #if !defined(HAVE_LIBPROJ) || !defined(HAVE_PROJECTS_H) 00504 00505 int GTIFProj4ToLatLong( GTIFDefn * psDefn, int nPoints, 00506 double *padfX, double *padfY ) 00507 { 00508 #ifdef DEBUG 00509 fprintf( stderr, 00510 "GTIFProj4ToLatLong() - PROJ.4 support not compiled in.\n" ); 00511 #endif 00512 return FALSE; 00513 } 00514 00515 int GTIFProj4FromLatLong( GTIFDefn * psDefn, int nPoints, 00516 double *padfX, double *padfY ) 00517 { 00518 #ifdef DEBUG 00519 fprintf( stderr, 00520 "GTIFProj4FromLatLong() - PROJ.4 support not compiled in.\n" ); 00521 #endif 00522 return FALSE; 00523 } 00524 #else 00525 00526 #include "projects.h" 00527 00528 #ifdef USE_PROJUV 00529 # define UV projUV 00530 #endif 00531 00532 /************************************************************************/ 00533 /* GTIFProj4FromLatLong() */ 00534 /* */ 00535 /* Convert lat/long values to projected coordinate for a */ 00536 /* particular definition. */ 00537 /************************************************************************/ 00538 00539 int GTIFProj4FromLatLong( GTIFDefn * psDefn, int nPoints, 00540 double *padfX, double *padfY ) 00541 00542 { 00543 char *pszProjection, **papszArgs; 00544 PJ *psPJ; 00545 int i; 00546 00547 /* -------------------------------------------------------------------- */ 00548 /* Get a projection definition. */ 00549 /* -------------------------------------------------------------------- */ 00550 pszProjection = GTIFGetProj4Defn( psDefn ); 00551 00552 if( pszProjection == NULL ) 00553 return FALSE; 00554 00555 /* -------------------------------------------------------------------- */ 00556 /* Parse into tokens for pj_init(), and initialize the projection. */ 00557 /* -------------------------------------------------------------------- */ 00558 00559 papszArgs = CSLTokenizeStringComplex( pszProjection, " +", TRUE, FALSE ); 00560 free( pszProjection ); 00561 00562 psPJ = pj_init( CSLCount(papszArgs), papszArgs ); 00563 00564 CSLDestroy( papszArgs ); 00565 00566 if( psPJ == NULL ) 00567 { 00568 return FALSE; 00569 } 00570 00571 /* -------------------------------------------------------------------- */ 00572 /* Process each of the points. */ 00573 /* -------------------------------------------------------------------- */ 00574 for( i = 0; i < nPoints; i++ ) 00575 { 00576 UV sUV; 00577 00578 sUV.u = padfX[i] * DEG_TO_RAD; 00579 sUV.v = padfY[i] * DEG_TO_RAD; 00580 00581 sUV = pj_fwd( sUV, psPJ ); 00582 00583 padfX[i] = sUV.u; 00584 padfY[i] = sUV.v; 00585 } 00586 00587 return TRUE; 00588 } 00589 00590 /************************************************************************/ 00591 /* GTIFProj4ToLatLong() */ 00592 /* */ 00593 /* Convert projection coordinates to lat/long for a particular */ 00594 /* definition. */ 00595 /************************************************************************/ 00596 00597 int GTIFProj4ToLatLong( GTIFDefn * psDefn, int nPoints, 00598 double *padfX, double *padfY ) 00599 00600 { 00601 char *pszProjection, **papszArgs; 00602 PJ *psPJ; 00603 int i; 00604 00605 /* -------------------------------------------------------------------- */ 00606 /* Get a projection definition. */ 00607 /* -------------------------------------------------------------------- */ 00608 pszProjection = GTIFGetProj4Defn( psDefn ); 00609 00610 if( pszProjection == NULL ) 00611 return FALSE; 00612 00613 /* -------------------------------------------------------------------- */ 00614 /* Parse into tokens for pj_init(), and initialize the projection. */ 00615 /* -------------------------------------------------------------------- */ 00616 00617 papszArgs = CSLTokenizeStringComplex( pszProjection, " +", TRUE, FALSE ); 00618 free( pszProjection ); 00619 00620 psPJ = pj_init( CSLCount(papszArgs), papszArgs ); 00621 00622 CSLDestroy( papszArgs ); 00623 00624 if( psPJ == NULL ) 00625 { 00626 return FALSE; 00627 } 00628 00629 /* -------------------------------------------------------------------- */ 00630 /* Process each of the points. */ 00631 /* -------------------------------------------------------------------- */ 00632 for( i = 0; i < nPoints; i++ ) 00633 { 00634 UV sUV; 00635 00636 sUV.u = padfX[i]; 00637 sUV.v = padfY[i]; 00638 00639 sUV = pj_inv( sUV, psPJ ); 00640 00641 padfX[i] = sUV.u * RAD_TO_DEG; 00642 padfY[i] = sUV.v * RAD_TO_DEG; 00643 } 00644 00645 return TRUE; 00646 } 00647 00648 00649 #endif /* has projects.h and -lproj */ 00650
1.2.3-20001105 written by Dimitri van Heesch,
© 1997-2000