[gdal-dev] how to add the gcp to a tif and make it to a geotiff?

[chao ding]

how to add the gcp to a tif and make it to a geotiff?
use gdal_translate -gcp we could do this.
(in gdal_translate.cpp=A3=A9
else if( EQUAL(argv[i],"-gcp") && i < argc - 4 )
{
   /* -gcp pixel line easting northing [elev] */

   nGCPCount++;
   pasGCPs =3D (GDAL_GCP *)
       CPLRealloc( pasGCPs, sizeof(GDAL_GCP) * nGCPCount );
   GDALInitGCPs( 1, pasGCPs + nGCPCount - 1 );


   pasGCPs[nGCPCount-1].dfGCPPixel =3D atof(argv[++i]);
   pasGCPs[nGCPCount-1].dfGCPLine =3D atof(argv[++i]);
   pasGCPs[nGCPCount-1].dfGCPX =3D atof(argv[++i]);
   pasGCPs[nGCPCount-1].dfGCPY =3D atof(argv[++i]);
   if( argv[i+1] !=3D NULL
       && (atof(argv[i+1]) !=3D 0.0 || argv[i+1][0] =3D=3D '0') )
       pasGCPs[nGCPCount-1].dfGCPZ =3D atof(argv[++i]);


   /* should set id and info? */



}


I wrtie   in my code =A3=BA

char* szFilter =3D "TIFF(*.tif)|*.tif|BMP (*.bmp)|*.bmp|JPG (*.jpg)|
*.jpg|Image(*.img)|*.img||";
CFileDialog     dlg(FALSE,TEXT(szFilter),NULL,OFN_HIDEREADONLY |
OFN_OVERWRITEPROMPT,TEXT(szFilter),this);

if(dlg.DoModal()=3D=3DIDOK)
{
       CString strSaveName =3D dlg.GetPathName();           //=BB=F1=B5=C3=
=B1=A3=B4=E6=B5=C4=CD=EA=D5=FB=C2=B7=BE=B6
       GDALAllRegister();
       GDALDataset* poDataset;
       poDataset =3D (GDALDataset *) GDALOpen(filename , GA_Update );

       if( !poDataset =3D=3D NULL )
       {
               VRTDataset *poVDS;//=B4=B4=BD=A8=D0=E9=C4=E2=CA=FD=BE=DD=BC=
=AF
               int nRasterXSize =3D poDataset->GetRasterXSize();
               int nRasterYSize =3D poDataset->GetRasterYSize();
               poVDS =3D (VRTDataset *) VRTCreate( nRasterXSize, nRasterYSi=
ze );
               char * pszOutputSRS =3D "WGS_84";
               poVDS->SetProjection( pszOutputSRS );//SRS=C9=E8=CE=AA=D0=C2=
=B5=BC=C8=EB=B5=C4SRS

               if(m_bChangeCoords)
               {
                       double        adfGeoTransform[6];
                       adfGeoTransform[0] =3D 20;        //m_dTLX;
                       adfGeoTransform[1] =3D 2; //m_dCellX;
                       adfGeoTransform[2] =3D 0;
                       adfGeoTransform[3] =3D 20;        //m_dTLY;
                       adfGeoTransform[4] =3D 0;
                       adfGeoTransform[5] =3D -2;        //m_dCellY;
                       poVDS->SetGeoTransform( adfGeoTransform );//        =
    }

               GDAL_GCP  * pasGCPs =3D NULL;
               int mGcpCount =3D 4;
               //gcplist
               pasGCPs =3D (GDAL_GCP *)CPLRealloc( pasGCPs, sizeof(GDAL_GCP=
) *
mGcpCount );

               pasGCPs[0].dfGCPPixel =3D 0;
               pasGCPs[0].dfGCPLine =3D 0;
               pasGCPs[0].dfGCPX =3D 863909;
               pasGCPs[0].dfGCPY =3D 993283;
               pasGCPs[0].dfGCPZ =3D 0;
               pasGCPs[0].pszId =3D "1";
               pasGCPs[0].pszInfo =3D "info 1";


               pasGCPs[1].dfGCPPixel =3D 0;
               pasGCPs[1].dfGCPLine =3D 201;
               pasGCPs[1].dfGCPX =3D 863909;
               pasGCPs[1].dfGCPY =3D 1024899;
               pasGCPs[1].dfGCPZ =3D 0;
               pasGCPs[1].pszId =3D "2";
               pasGCPs[1].pszInfo =3D "info 2";

               pasGCPs[2].dfGCPPixel =3D 201;
               pasGCPs[2].dfGCPLine =3D 201;
               pasGCPs[2].dfGCPX =3D 809988;
               pasGCPs[2].dfGCPY =3D 1024899;
               pasGCPs[2].dfGCPZ =3D 0;
               pasGCPs[2].pszId =3D "3";
               pasGCPs[2].pszInfo =3D "info 3";

               pasGCPs[3].dfGCPPixel =3D 201;
               pasGCPs[3].dfGCPLine =3D 0;
               pasGCPs[3].dfGCPX =3D 809988;
               pasGCPs[3].dfGCPY =3D 0;
               pasGCPs[3].dfGCPZ =3D 0;
               pasGCPs[3].pszId =3D "4";
               pasGCPs[3].pszInfo =3D "info 4";

               char * m_GcpProjection =3D "WGS_84";
               poVDS->SetGCPs(mGcpCount,pasGCPs,m_GcpProjection);
               if
(poVDS->SetGCPs(mGcpCount,pasGCPs,m_GcpProjection)=3D=3DCE_None) //=CB=B5=
=C3=F7
setgcp=CA=C7=B3=C9=B9=A6=B5=C4,=B5=AB=CA=C7=D5=E2=C0=EF=B2=A2=CE=B4=B8=FC=
=B8=C4tif=A3=AC
               {
                       AfxMessageBox("ok");
               }

               double madfGeoTransform[6] =3D {444720, 2.5, 0, 3751320, 0, =
-2.5};/
               char * src_adfGeoTransform =3D poDataset->GetGeoTransform();
               poVDS->SetGeoTransform( madfGeoTransform);
               poVDS->SetProjection(mProjection);
               poVDS->SetGCPs(mGcpCount,pasGCPs,m_GcpProjection);
               poVDS->SetMetadata( poDataset->GetMetadata( ) );//=D4=AA=CA=
=FD=BE=DD=B1=A3=B3=D6=BE=C9=B5=C4

               char * mGcpProjection;
               poVDS->GetGCPProjection();

               AfxMessageBox(poVDS->GetGCPProjection());
               AfxMessageBox(poVDS->GetGCPCount());
               AfxMessageBox(m_GcpProjection);

               int i =3D 0;
               int nBandCount =3D poDataset->GetRasterCount();
               for( i =3D 0; i < nBandCount; i++ )
               {
                       VRTSourcedRasterBand   *poVRTBand;
                       GDALRasterBand  *poSrcBand;
                       GDALDataType    eBandType;

                       poSrcBand =3D poDataset->GetRasterBand(i+1);
                       eBandType =3D poSrcBand->GetRasterDataType();
                       //       Create this band.
                       poVDS->AddBand( eBandType, NULL );
                       poVRTBand =3D (VRTSourcedRasterBand *)
poVDS->GetRasterBand( i+1 );
                       poVRTBand->AddSimpleSource( poSrcBand, 0, 0,nRasterX=
Size,
nRasterYSize,0,0,nRasterXSize, nRasterYSize);


                       /*
--------------------------------------------------------------------
*/
                       /*      copy over some other information of
interest.    */
                       /*
--------------------------------------------------------------------
*/
                       poVRTBand->CopyCommonInfoFrom( poSrcBand );
               }//for
               //  Write to the output file using CopyCreate().
               GDALDriver * poDriver =3D NULL;
               poDriver =3D poDataset->GetDriver();
               GDALDataset* pNewDs =3D poDriver->CreateCopy( (LPCTSTR)strSa=
veName,
poVDS, FALSE, NULL,NULL,NULL);
               //GDALDataset* mNewDs =3D RasterIO()
               if( pNewDs!=3D NULL)
                       delete pNewDs;

               delete poVDS;
               CPLFree( pszOutputSRS );


       }//if !poDataset =3D=3D NULL

       if(poDataset !=3D NULL)
       delete poDataset;