erdas img文件如何取高程呢？

[xiongyouyi]

我现在想用GDAL/OGR从shapefile中读取point的xy坐标，然后从img中找z，
对img格式不太了解，请问已知xy如何取z呢？谢谢
找了一段代码，不知道这个是什么意思，
p3dDEMGridPoint[index][2] = pScanline[x]/200;//pScanline[x]+offfset;

#include "HFADataset.h"     
bool ImportDEMDataERDASImg::LoadSingleFile(LPCSTR lpszPathName)
{
         HFADataset* pDataset = Open( lpszPathName );  //通过文件路径打开pDataSet,来源于GDAL
         if( pDataset != NULL )
        {
               pDataset->GetProjectionRef();

              double        adfGeoTransform[6];
              pDataset->GetGeoTransform( adfGeoTransform );

              double Cellsize;           //单元格大小
             double xMin = 0.0;       //X轴最小数值
             double yMin = 0.0;       //Y轴最小数值

              int  x_pos, y_pos;
              double GeoX, GeoY; 
              int  x, y;
              long index = 0;
              for(int i=0 ; i< pDataset->GetRasterCount(); i++ )
              {
                   if(adfGeoTransform[5]<0.0)                //坐标变换
                        GdalToWorld(adfGeoTransform, 0 , pDataset->GetRasterYSize(), GeoX,  GeoY );
                   else
                   {
                        adfGeoTransform[5]*=-1.0;
                        GdalToWorld(adfGeoTransform, 0 , pDataset->GetRasterYSize(), GeoX,  GeoY );
                   }

                   Cellsize = adfGeoTransform[1];
                   xMin = GeoX;
                   yMin = GeoY;

                   HFARasterBand  *pBand;
                   pBand = pDataset->GetRasterBand( i+1 );          //从1开始计数，具体的波段

                  double offfset = pBand->GetOffset();          //波段中的偏移量

                  float * pScanline = (float *) malloc(sizeof(float)*pDataset->GetRasterXSize());  //沿X轴长度的扫描线，
 
                  int iRows =  pDataset->GetRasterYSize();    //行数
                  int iColumns = pDataset->GetRasterXSize();//列数

                   double dMin = 99999;
                   double dMax = -99999;               //准备记录高程数值的最小和最大数值
                   POINT3d* p3dDEMGridPoint = new POINT3d[iRows*iColumns];      //added by houtao 2007 10 29
                   for (y = 0; y < iRows; y++)
                   {

                        pBand->RasterIO( GF_Read, 0, y, pDataset->GetRasterXSize(), 1, 
                                                    pScanline, pDataset->GetRasterXSize(), 1, GDT_Float32, 
                                                     0, 0 );                //读取数据到pScanline中

                    for (x = 0; x < iColumns; x++)
                    { 
                         GdalToWorld(adfGeoTransform, x , y, GeoX,  GeoY );
                         WorldToGeoView (xMin, yMin,Cellsize, GeoX, GeoY,x_pos, y_pos );
     
                          index = y*iColumns + x;
                         p3dDEMGridPoint[index][0] = x_pos;
                         p3dDEMGridPoint[index][1] = y_pos;                //记录数据
                         p3dDEMGridPoint[index][2] = pScanline[x]/200;//pScanline[x]+offfset;

                         if( p3dDEMGridPoint[index][2] > -500 && p3dDEMGridPoint[index][2] < 9000) //高程数值在正常的范围内
                         {                                            //寻找最小数值
                              if( dMin > p3dDEMGridPoint[index][2] ) dMin = p3dDEMGridPoint[index][2];
                              if( dMax < p3dDEMGridPoint[index][2] ) dMax = p3dDEMGridPoint[index][2]; 
                         }
                       }                //对X轴循环结束
                   }                    //对Y轴循环结束

                   m_p3dDEMGrid->SetRows(iRows);          //记下行数
                   m_p3dDEMGrid->SetColumns(iColumns);         //记下列数
                   m_p3dDEMGrid->SetSurfacePoint( p3dDEMGridPoint );     //设置点集
                   m_p3dDEMGrid->SetZMin(dMin);          //记下Z轴最小数值
                   m_p3dDEMGrid->SetZMax(dMax);          //记下Z轴最大数值
                   m_p3dDEMGrid->SetXMin(xMin);          //记下X轴最小数值
                   m_p3dDEMGrid->SetXMax(xMin+iColumns*Cellsize);      //记下X轴最大数值
                   m_p3dDEMGrid->SetYMin(yMin);          //记下Y轴最小数值
                   m_p3dDEMGrid->SetYMax(yMin+iRows*Cellsize);       //记下Y轴最大数值
                   free(pScanline);
            }                   //结束这个波段的循环

         }

         return true;
}