hi adam,
should be an easy fix in that you can do in your GDAL code look how
the projection setting is done in GeoProjectionConverter.cpp from
LASlib. There it only sets up a user defined projection if all
parameters for that projection were provided. This happens in the
set_projection_from_geo_keys() function (see below).
cheers,
martin
bool GeoProjectionConverter::set_projection_from_geo_keys(int
num_geo_keys, GeoProjectionGeoKeys* geo_keys, char* geo_ascii_params,
double* geo_double_params, char* description)
{
bool user_defined_ellipsoid = false;
int user_defined_projection = 0;
int offsetProjStdParallel1GeoKey = -1;
int offsetProjStdParallel2GeoKey = -1;
int offsetProjNatOriginLatGeoKey = -1;
int offsetProjFalseEastingGeoKey = -1;
int offsetProjFalseNorthingGeoKey = -1;
int offsetProjCenterLongGeoKey = -1;
int offsetProjScaleAtNatOriginGeoKey = -1;
bool has_projection = false;
int ellipsoid = -1;
this->num_geo_keys = num_geo_keys;
this->geo_keys = geo_keys;
this->geo_ascii_params = geo_ascii_params;
this->geo_double_params = geo_double_params;
for (int i = 0; i < num_geo_keys; i++)
{
switch (geo_keys[i].key_id)
{
case 1024: // GTModelTypeGeoKey
if (geo_keys[i].value_offset == 2) // ModelTypeGeographic
{
has_projection = set_longlat_projection(description);
}
break;
case 2048: // GeographicTypeGeoKey
switch (geo_keys[i].value_offset)
{
case 32767: // user-defined GCS
user_defined_ellipsoid = true;
break;
case 4001: // GCSE_Airy1830
ellipsoid = 1;
break;
case 4002: // GCSE_AiryModified1849
ellipsoid = 16;
break;
case 4003: // GCSE_AustralianNationalSpheroid
ellipsoid = 2;
break;
case 4004: // GCSE_Bessel1841
case 4005: // GCSE_Bessel1841Modified
ellipsoid = 3;
break;
case 4006: // GCSE_BesselNamibia
ellipsoid = 4;
break;
case 4008: // GCSE_Clarke1866
case 4009: // GCSE_Clarke1866Michigan
ellipsoid = GEO_ELLIPSOID_NAD27;
break;
case 4010: // GCSE_Clarke1880_Benoit
case 4011: // GCSE_Clarke1880_IGN
case 4012: // GCSE_Clarke1880_RGS
case 4013: // GCSE_Clarke1880_Arc
case 4014: // GCSE_Clarke1880_SGA1922
case 4034: // GCSE_Clarke1880
ellipsoid = 6;
break;
case 4015: // GCSE_Everest1830_1937Adjustment
case 4016: // GCSE_Everest1830_1967Definition
case 4017: // GCSE_Everest1830_1975Definition
ellipsoid = 7;
break;
case 4018: // GCSE_Everest1830Modified
ellipsoid = 17;
break;
case 4019: // GCSE_GRS1980
ellipsoid = GEO_ELLIPSOID_NAD83;
break;
case 4020: // GCSE_Helmert1906
ellipsoid = 12;
break;
case 4022: // GCSE_International1924
case 4023: // GCSE_International1967
ellipsoid = 14;
break;
case 4024: // GCSE_Krassowsky1940
ellipsoid = 15;
break;
case 4030: // GCSE_WGS84
ellipsoid = GEO_ELLIPSOID_WGS84;
break;
case 4267: // GCS_NAD27
ellipsoid = GEO_ELLIPSOID_NAD27;
break;
case 4269: // GCS_NAD83
ellipsoid = GEO_ELLIPSOID_NAD83;
break;
case 4322: // GCS_WGS_72
ellipsoid = GEO_ELLIPSOID_WGS72;
break;
case 4326: // GCS_WGS_84
ellipsoid = GEO_ELLIPSOID_WGS84;
break;
default:
fprintf(stderr, "GeographicTypeGeoKey: look-up for %d not
implemented\n", geo_keys[i].value_offset);
}
break;
case 2050: // GeogGeodeticDatumGeoKey
switch (geo_keys[i].value_offset)
{
case 32767: // user-defined GCS
user_defined_ellipsoid = true;
break;
case 6202: // Datum_Australian_Geodetic_Datum_1966
case 6203: // Datum_Australian_Geodetic_Datum_1984
ellipsoid = 2;
break;
case 6267: // Datum_North_American_Datum_1927
ellipsoid = GEO_ELLIPSOID_NAD27;
break;
case 6269: // Datum_North_American_Datum_1983
ellipsoid = GEO_ELLIPSOID_NAD83;
break;
case 6322: // Datum_WGS72
ellipsoid = GEO_ELLIPSOID_WGS72;
break;
case 6326: // Datum_WGS84
ellipsoid = GEO_ELLIPSOID_WGS84;
break;
case 6001: // DatumE_Airy1830
ellipsoid = 1;
break;
case 6002: // DatumE_AiryModified1849
ellipsoid = 16;
break;
case 6003: // DatumE_AustralianNationalSpheroid
ellipsoid = 2;
break;
case 6004: // DatumE_Bessel1841
case 6005: // DatumE_BesselModified
ellipsoid = 3;
break;
case 6006: // DatumE_BesselNamibia
ellipsoid = 4;
break;
case 6008: // DatumE_Clarke1866
case 6009: // DatumE_Clarke1866Michigan
ellipsoid = GEO_ELLIPSOID_NAD27;
break;
case 6010: // DatumE_Clarke1880_Benoit
case 6011: // DatumE_Clarke1880_IGN
case 6012: // DatumE_Clarke1880_RGS
case 6013: // DatumE_Clarke1880_Arc
case 6014: // DatumE_Clarke1880_SGA1922
case 6034: // DatumE_Clarke1880
ellipsoid = 6;
break;
case 6015: // DatumE_Everest1830_1937Adjustment
case 6016: // DatumE_Everest1830_1967Definition
case 6017: // DatumE_Everest1830_1975Definition
ellipsoid = 7;
break;
case 6018: // DatumE_Everest1830Modified
ellipsoid = 17;
break;
case 6019: // DatumE_GRS1980
ellipsoid = GEO_ELLIPSOID_NAD83;
break;
case 6020: // DatumE_Helmert1906
ellipsoid = 12;
break;
case 6022: // DatumE_International1924
case 6023: // DatumE_International1967
ellipsoid = 14;
break;
case 6024: // DatumE_Krassowsky1940
ellipsoid = 15;
break;
case 6030: // DatumE_WGS84
ellipsoid = GEO_ELLIPSOID_WGS84;
break;
default:
fprintf(stderr, "GeogGeodeticDatumGeoKey: look-up for %d not
implemented\n", geo_keys[i].value_offset);
}
break;
case 2052: // GeogLinearUnitsGeoKey
switch (geo_keys[i].value_offset)
{
case 9001: // Linear_Meter
set_coordinates_in_meter();
break;
case 9002: // Linear_Foot
set_coordinates_in_feet();
break;
case 9003: // Linear_Foot_US_Survey
set_coordinates_in_survey_feet();
break;
default:
fprintf(stderr, "GeogLinearUnitsGeoKey: look-up for %d not
implemented\n", geo_keys[i].value_offset);
}
break;
case 2056: // GeogEllipsoidGeoKey
switch (geo_keys[i].value_offset)
{
case 7001: // Ellipse_Airy_1830
ellipsoid = 1;
break;
case 7002: // Ellipse_Airy_Modified_1849
ellipsoid = 16;
break;
case 7003: // Ellipse_Australian_National_Spheroid
ellipsoid = 2;
break;
case 7004: // Ellipse_Bessel_1841
case 7005: // Ellipse_Bessel_Modified
ellipsoid = 3;
break;
case 7006: // Ellipse_Bessel_Namibia
ellipsoid = 4;
break;
case 7008: // Ellipse_Clarke_1866
case 7009: // Ellipse_Clarke_1866_Michigan
ellipsoid = GEO_ELLIPSOID_NAD27;
break;
case 7010: // Ellipse_Clarke1880_Benoit
case 7011: // Ellipse_Clarke1880_IGN
case 7012: // Ellipse_Clarke1880_RGS
case 7013: // Ellipse_Clarke1880_Arc
case 7014: // Ellipse_Clarke1880_SGA1922
case 7034: // Ellipse_Clarke1880
ellipsoid = 6;
break;
case 7015: // Ellipse_Everest1830_1937Adjustment
case 7016: // Ellipse_Everest1830_1967Definition
case 7017: // Ellipse_Everest1830_1975Definition
ellipsoid = 7;
break;
case 7018: // Ellipse_Everest1830Modified
ellipsoid = 17;
break;
case 7019: // Ellipse_GRS_1980
ellipsoid = GEO_ELLIPSOID_NAD83;
break;
case 7020: // Ellipse_Helmert1906
ellipsoid = 12;
break;
case 7022: // Ellipse_International1924
case 7023: // Ellipse_International1967
ellipsoid = 14;
break;
case 7024: // Ellipse_Krassowsky1940
ellipsoid = 15;
break;
case 7030: // Ellipse_WGS_84
ellipsoid = GEO_ELLIPSOID_WGS84;
break;
default:
fprintf(stderr, "GeogEllipsoidGeoKey: look-up for %d not
implemented\n", geo_keys[i].value_offset);
}
break;
case 3072: // ProjectedCSTypeGeoKey
if (geo_keys[i].value_offset != 32767)
has_projection =
set_ProjectedCSTypeGeoKey(geo_keys[i].value_offset, description);
break;
case 3075: // ProjCoordTransGeoKey
user_defined_projection = 0;
switch (geo_keys[i].value_offset)
{
case 1: // CT_TransverseMercator
user_defined_projection = 1;
break;
case 8: // CT_LambertConfConic_2SP
user_defined_projection = 8;
break;
case 2: // CT_TransvMercator_Modified_Alaska
fprintf(stderr, "ProjCoordTransGeoKey:
CT_TransvMercator_Modified_Alaska not implemented\n");
break;
case 3: // CT_ObliqueMercator
fprintf(stderr, "ProjCoordTransGeoKey: CT_ObliqueMercator not
implemented\n");
break;
case 4: // CT_ObliqueMercator_Laborde
fprintf(stderr, "ProjCoordTransGeoKey:
CT_ObliqueMercator_Laborde not implemented\n");
break;
case 5: // CT_ObliqueMercator_Rosenmund
fprintf(stderr, "ProjCoordTransGeoKey:
CT_ObliqueMercator_Rosenmund not implemented\n");
break;
case 6: // CT_ObliqueMercator_Spherical
fprintf(stderr, "ProjCoordTransGeoKey:
CT_ObliqueMercator_Spherical not implemented\n");
break;
case 7: // CT_Mercator
fprintf(stderr, "ProjCoordTransGeoKey: CT_Mercator not implemented\n");
break;
case 9: // CT_LambertConfConic_Helmert
fprintf(stderr, "ProjCoordTransGeoKey:
CT_LambertConfConic_Helmert not implemented\n");
break;
case 10: // CT_LambertAzimEqualArea
fprintf(stderr, "ProjCoordTransGeoKey: CT_LambertAzimEqualArea
not implemented\n");
break;
case 11: // CT_AlbersEqualArea
fprintf(stderr, "ProjCoordTransGeoKey: CT_AlbersEqualArea not
implemented\n");
break;
case 12: // CT_AzimuthalEquidistant
fprintf(stderr, "ProjCoordTransGeoKey: CT_AzimuthalEquidistant
not implemented\n");
break;
case 13: // CT_EquidistantConic
fprintf(stderr, "ProjCoordTransGeoKey: CT_EquidistantConic not
implemented\n");
break;
case 14: // CT_Stereographic
fprintf(stderr, "ProjCoordTransGeoKey: CT_Stereographic not
implemented\n");
break;
case 15: // CT_PolarStereographic
fprintf(stderr, "ProjCoordTransGeoKey: CT_PolarStereographic
not implemented\n");
break;
case 16: // CT_ObliqueStereographic
fprintf(stderr, "ProjCoordTransGeoKey: CT_ObliqueStereographic
not implemented\n");
break;
case 17: // CT_Equirectangular
fprintf(stderr, "ProjCoordTransGeoKey: CT_Equirectangular not
implemented\n");
break;
case 18: // CT_CassiniSoldner
fprintf(stderr, "ProjCoordTransGeoKey: CT_CassiniSoldner not
implemented\n");
break;
case 19: // CT_Gnomonic
fprintf(stderr, "ProjCoordTransGeoKey: CT_Gnomonic not implemented\n");
break;
case 20: // CT_MillerCylindrical
fprintf(stderr, "ProjCoordTransGeoKey: CT_MillerCylindrical
not implemented\n");
break;
case 21: // CT_Orthographic
fprintf(stderr, "ProjCoordTransGeoKey: CT_Orthographic not
implemented\n");
break;
case 22: // CT_Polyconic
fprintf(stderr, "ProjCoordTransGeoKey: CT_Polyconic not implemented\n");
break;
case 23: // CT_Robinson
fprintf(stderr, "ProjCoordTransGeoKey: CT_Robinson not implemented\n");
break;
case 24: // CT_Sinusoidal
fprintf(stderr, "ProjCoordTransGeoKey: CT_Sinusoidal not
implemented\n");
break;
case 25: // CT_VanDerGrinten
fprintf(stderr, "ProjCoordTransGeoKey: CT_VanDerGrinten not
implemented\n");
break;
case 26: // CT_NewZealandMapGrid
fprintf(stderr, "ProjCoordTransGeoKey: CT_NewZealandMapGrid
not implemented\n");
break;
case 27: // CT_TransvMercator_SouthOriented
fprintf(stderr, "ProjCoordTransGeoKey:
CT_TransvMercator_SouthOriented not implemented\n");
break;
default:
fprintf(stderr, "ProjCoordTransGeoKey: look-up for %d not
implemented\n", geo_keys[i].value_offset);
}
break;
case 3076: // ProjLinearUnitsGeoKey
set_ProjLinearUnitsGeoKey(geo_keys[i].value_offset);
break;
case 3078: // ProjStdParallel1GeoKey
offsetProjStdParallel1GeoKey = geo_keys[i].value_offset;
break;
case 3079: // ProjStdParallel2GeoKey
offsetProjStdParallel2GeoKey = geo_keys[i].value_offset;
break;
case 3081: // ProjNatOriginLatGeoKey
offsetProjNatOriginLatGeoKey = geo_keys[i].value_offset;
break;
case 3082: // ProjFalseEastingGeoKey
offsetProjFalseEastingGeoKey = geo_keys[i].value_offset;
break;
case 3083: // ProjFalseNorthingGeoKey
offsetProjFalseNorthingGeoKey = geo_keys[i].value_offset;
break;
case 3088: // ProjCenterLongGeoKey
offsetProjCenterLongGeoKey = geo_keys[i].value_offset;
break;
case 3092: // ProjScaleAtNatOriginGeoKey
offsetProjScaleAtNatOriginGeoKey = geo_keys[i].value_offset;
break;
case 4096: // VerticalCSTypeGeoKey
set_VerticalCSTypeGeoKey(geo_keys[i].value_offset);
break;
case 4099: // VerticalUnitsGeoKey
set_VerticalUnitsGeoKey(geo_keys[i].value_offset);
break;
}
}
if (ellipsoid != -1)
{
set_reference_ellipsoid(ellipsoid);
}
if (!has_projection)
{
if (user_defined_projection == 1)
{
if ((offsetProjFalseEastingGeoKey >= 0) &&
(offsetProjFalseNorthingGeoKey >= 0) &&
(offsetProjNatOriginLatGeoKey >= 0) &&
(offsetProjCenterLongGeoKey >= 0) &&
(offsetProjScaleAtNatOriginGeoKey >= 0))
{
double falseEastingMeter =
geo_double_params[offsetProjFalseEastingGeoKey] * coordinates2meter;
double falseNorthingMeter =
geo_double_params[offsetProjFalseNorthingGeoKey] * coordinates2meter;
double latOriginDeg = geo_double_params[offsetProjNatOriginLatGeoKey];
double longMeridianDeg = geo_double_params[offsetProjCenterLongGeoKey];
double scaleFactor =
geo_double_params[offsetProjScaleAtNatOriginGeoKey];
set_transverse_mercator_projection(falseEastingMeter,
falseNorthingMeter, latOriginDeg, longMeridianDeg, scaleFactor);
if (description)
{
sprintf(description, "generic traverse mercator");
}
has_projection = true;
}
}
else if (user_defined_projection == 8)
{
if ((offsetProjFalseEastingGeoKey >= 0) &&
(offsetProjFalseNorthingGeoKey >= 0) &&
(offsetProjNatOriginLatGeoKey >= 0) &&
(offsetProjCenterLongGeoKey >= 0) &&
(offsetProjStdParallel1GeoKey >= 0) &&
(offsetProjStdParallel2GeoKey >= 0))
{
double falseEastingMeter =
geo_double_params[offsetProjFalseEastingGeoKey] * coordinates2meter;
double falseNorthingMeter =
geo_double_params[offsetProjFalseNorthingGeoKey] * coordinates2meter;
double latOriginDeg = geo_double_params[offsetProjNatOriginLatGeoKey];
double longOriginDeg = geo_double_params[offsetProjCenterLongGeoKey];
double firstStdParallelDeg =
geo_double_params[offsetProjStdParallel1GeoKey];
double secondStdParallelDeg =
geo_double_params[offsetProjStdParallel2GeoKey];
set_lambert_conformal_conic_projection(falseEastingMeter,
falseNorthingMeter, latOriginDeg, longOriginDeg, firstStdParallelDeg,
secondStdParallelDeg);
if (description)
{
sprintf(description, "generic lambert conformal conic");
}
has_projection = true;
}
}
}
return has_projection;
}
_______________________________________________
Liblas-devel mailing list
[email protected]
http://lists.osgeo.org/mailman/listinfo/liblas-devel