Take the 2-minute tour ×
Geographic Information Systems Stack Exchange is a question and answer site for cartographers, geographers and GIS professionals. It's 100% free, no registration required.

I'm using libgeotiff and proj4 to read and (try to!) process GeoTIFFs. After reading the file, I want to get the corner coordinates in WGS84/UTM - but there are so many functions and parameters to try and understand, I'm flailing around a bit. A lot actually!

My current hacky function goes like this:

//------------------------------------------------------------------------------
BOOL CGeoTIFFImporter::GetUTMCorner(Corners cnr, double& x, double& y)
//------------------------------------------------------------------------------
{
    BOOL bRet = FALSE;
    if (m_bIsGeoTIFF)
    {
        x = -DBL_MAX;
        y = -DBL_MAX;
        double xtmp = -DBL_MAX;
        double ytmp = -DBL_MAX;
        CString spjUTM = "";
        projPJ pjUTM, pjPJ;
        char *pszProjection;

        switch (cnr)
        {
            case LowerLeft: x = 0.0; y = m_Ypx; break;
            case LowerRight: x = m_Xpx; y = m_Ypx; break;
            case UpperLeft: x = 0.0; y = 0.0; break;
            case UpperRight: x = m_Xpx; y = 0.0; break;
        }

        // Try to transform the coordinate into PCS space
        bRet = GTIFImageToPCS(m_pGTF, &x, &y);

        if (bRet && m_GTFDef.Zone != 0)
        {
            xtmp = x; ytmp = y;
            if (m_GTFDef.Model == ModelTypeGeographic)
            {
            }
            else if (!GTIFProj4ToLatLong(&m_GTFDef, 1, &xtmp, &ytmp))
            {
                return FALSE;
            }

            spjUTM.Format("+proj=utm +zone=%d +ellps=WGS84 +datum=WGS84", UTMZone(xtmp));

            pjUTM = pj_init_plus(spjUTM);

            if (pjUTM == NULL)
            {
                return FALSE;
            }

            pszProjection = GTIFGetProj4Defn(&m_GTFDef);

            if (pszProjection == NULL)
            {
                return FALSE;
            }
            pjPJ = pj_init_plus(pszProjection);
            free(pszProjection);

            if (pjPJ == NULL)
            {
                return FALSE;
            }

            bRet = (pj_transform(pjPJ, pjUTM, 1, 1, &x, &y, NULL) == 0);
        }
    }
    return bRet;
}

Somewhat unsurprisingly given my state of knowledge onthe subject, this doesn't work. Can anyone give me some pointers as to what I'm doing wrong?

Cheers.

share|improve this question
    
Typical! I've been Googling this for a week without success - as soon as I post this, I find the answer - I think! Just testing (ytmp > 0) and adding +north or +south as appropriate. Still, I'll leave this open for comments for a while, because I'm not 100% sure! –  Kyudos Jun 18 '13 at 22:20
    
What are you actually trying to do? Hopefully not a transformation, as this can easily be done in GDAL: gdal.org/gdalwarp.html (as I'm sure you're aware!) –  alexgleith Jun 18 '13 at 23:43
    
I'm importing and placing GeoTIFFs in my CAD app, but I want to stick to UTM and WGS84 for simplicity and planar-ness. GDAL is great, but it has a massive overhead for this one function I actually want. –  Kyudos Jun 18 '13 at 23:49
1  
So you are writing an app and are writing GeoTIFF reading functionality? Fair enough. But many applications use GDAL for this and it means that you can load more than just GeoTIFFs: trac.osgeo.org/gdal/wiki/SoftwareUsingGdal –  alexgleith Jun 19 '13 at 0:13
2  
That's true, but several tens of MBs of DLLs to support functions I'm don't (currently) want or need is just too bloaty for me! If I can add GeoTIFF and KML support, that'll do me for now. –  Kyudos Jun 19 '13 at 0:30
add comment

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.