C++ GDALCreate函数代码示例

gint gv_view_area_print_to_file(GvViewArea *view, int width, int height, const char * filename, const char * format, int is_rgb)
{
    GDALDriverH   driver;
    GDALDatasetH  dataset;
    gint          errcode;

    driver = GDALGetDriverByName( format );
    if( driver == NULL )
        return -1;

    if( is_rgb )
        dataset = GDALCreate( driver, filename, width, height, 3, GDT_Byte, 
                              NULL );
    else
        dataset = GDALCreate( driver, filename, width, height, 1, GDT_Byte, 
                              NULL );

    if( dataset == NULL )
        return -1;

    errcode = gv_view_area_render_to_func( view, width, height, 
                                           print_handler, dataset );
    GDALClose( dataset );

    print_handler( NULL, NULL );
    
    return errcode;
}

/* Makes a copy of a dataset, and opens it for writing.. */
GDALDatasetH
make_me_a_sandwitch (GDALDatasetH * in_dataset, const char *copy_file_name)
{
  char **papszOptions = NULL;
  const char *pszFormat = "GTiff";
  double adfGeoTransform[6];
  GDALDriverH hDriver;
  GDALDatasetH out_gdalfile;
  hDriver = GDALGetDriverByName (pszFormat);
  papszOptions = CSLSetNameValue (papszOptions, "TILED", "YES");
  papszOptions = CSLSetNameValue (papszOptions, "COMPRESS", "DEFLATE");

  /*Perhaps controversal - default to bigtiff... */
  /*papszOptions = CSLSetNameValue( papszOptions, "BIGTIFF", "YES" ); */

  /*return GDALCreateCopy( hDriver, copy_file_name, *in_dataset, FALSE, papszOptions, NULL, NULL ); */
  out_gdalfile = GDALCreate (hDriver, copy_file_name,
			     GDALGetRasterXSize (*in_dataset),
			     GDALGetRasterYSize (*in_dataset),
			     GDALGetRasterCount (*in_dataset),
			     GDT_Byte, papszOptions);

  /* Set geotransform */
  GDALGetGeoTransform (*in_dataset, adfGeoTransform);
  GDALSetGeoTransform (out_gdalfile, adfGeoTransform);

  /* Set projection */
  GDALSetProjection (out_gdalfile, GDALGetProjectionRef (*in_dataset));
  return out_gdalfile;
}

int sat_save_ch(s_sat * sat, const char * fname, unsigned ch)
{
	try;

		int ret = 0;
		unsigned height = sat->height, width = sat->width;
		GDALDatasetH ds = NULL;

		throw_null((ds = GDALCreate(drv_tiff, fname, width, height, 1, GDT_Byte, NULL)));

		throw((GDALSetProjection(ds, sat->proj_ref) == CE_Failure));
		throw((GDALSetGeoTransform(ds, sat->gt_coef) == CE_Failure));

		throw(GDALDatasetRasterIO(ds, GF_Write, 0, 0, width, height, sat->pixel[ch], width, height, GDT_Byte, 1, NULL, 0, 0, 0) == CE_Failure);

	catch;

		ret = -1;

	finally;

		if(ds != NULL)
			GDALClose(ds);

	return ret;
}

void write_image_file(const char * filename, int nx, int ny, float * image)
{
    int i;
    int nBands = 1;
    printf("writing image file %s (%d,%d)\n", filename, nx, ny);

    unsigned char * buf = (unsigned char*) malloc(nx*ny*sizeof(unsigned char));

    for (i = 0; i < nx*ny; i++) {
        buf[i] = (unsigned char) image[i];
    }

    GDALAllRegister();

    GDALDriverH driver = GDALGetDriverByName("GTiff");
    if (driver == NULL) {
        exit(1);
    }

    GDALDatasetH dataset = GDALCreate(driver, filename, nx, ny, nBands, GDT_Byte, NULL);
    if (dataset == NULL) {
        fprintf(stderr, "write_image_file: failed to open file %s\n", filename);
    }

    GDALDatasetRasterIO(dataset, GF_Write, 0, 0, nx, ny, buf, nx, ny, GDT_Byte,
                        nBands, NULL, 1, nx, nx*ny);

    GDALClose(dataset);
    free(buf);
}

//void CUtils::createNewByteGeoTIFF(const char* fileName, int bands, int rows, int cols, double adfGeoTransform[6], char szProjection[512], byte fillData, byte noDataValue)
void CUtils::createNewByteGeoTIFF(const char* fileName, int bands, int rows, int cols, double adfGeoTransform[6], const char * szProjection, byte fillData, byte noDataValue)
{
	char **papszOptions = NULL;
    GDALDriverH hDriver;
    GDALRasterBandH hDataset;
//	GDALRasterBandH hBand;
	if( (hDriver = GDALGetDriverByName("GTiff")) != NULL)
	{
		fprintf(stderr, "Create image %s...\n", fileName);
        if( (hDataset = GDALCreate( hDriver, fileName, cols, rows, bands, GDT_Byte, papszOptions )) !=NULL )
		{
			GDALSetGeoTransform(hDataset, adfGeoTransform );
			GDALSetProjection(hDataset, szProjection );
			/*
			int pr = CUtils::progress_ln_ex(stderr, 0, 0, START_PROGRESS);
			for(int band = 1; band<=bands; band++)
			{
				if( (hBand = GDALGetRasterBand(hDataset, band)) != NULL )
				{
					byte *pline = (byte *)CPLMalloc(sizeof(byte)*cols);
					for(int i=0; i<cols; i++) pline[i] = fillData;

					for(int i=0; i<rows; i++) GDALRasterIO(hBand, GF_Write, 0, i, cols, 1, pline, cols, 1, GDT_Byte, 0, 0 );
					CPLFree(pline);

					GDALSetRasterNoDataValue(hBand, noDataValue);
				}
				pr = CUtils::progress_ln_ex(stderr, band-1, bands, pr);
			}
			CUtils::progress_ln_ex(stderr, 0, 0, END_PROGRESS);
			*/
			GDALClose(hDataset);
		}
	}
}

//辐射校正处理=====================================================================================================================================================================================
//绝对辐射校正
long QPDLevel1Process::Level1Proc_RadiationAbsolute(const char* pathImg, const char* pathImgRad, const char* pathAbsRegFile)
{
	CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO");	//中文路径
	GDALAllRegister();
	long lError = 0;
	unsigned short *imgBuffer = NULL;	//影像数据
	float* parametersA = NULL, *parametersB = NULL, *parametersAux = NULL;//校正系数

	GDALDatasetH m_dataset = GDALOpen(pathImg, GA_ReadOnly);
	int xsize = GDALGetRasterXSize(m_dataset);
	int ysize = GDALGetRasterYSize(m_dataset);
	int bands = GDALGetRasterCount(m_dataset);

	char **papszOptions = NULL;
	papszOptions = CSLSetNameValue(papszOptions, "INTERLEAVE", "BAND");
	GDALDatasetH m_datasetdst = GDALCreate(GDALGetDriverByName("GTiff"), pathImgRad, xsize, ysize, bands, GDT_UInt16, papszOptions);

	//int nSamples, nLines, nLevels;
	//LevelProc_GetParameterInfo(pathImgRad, nSamples, nLines, nLevels);

	try
	{
		parametersA = new float[bands];
		parametersB = new float[bands];
		imgBuffer = new unsigned short[xsize*ysize];
	}
	catch (bad_alloc)
	{
		printf("allocate memory error\n");
		exit(-1);
	}

	Level1Proc_AbsoluteParameters(pathAbsRegFile, parametersA, parametersB);
	for (int i = 0; i < bands; i++)
	{
		GDALRasterIO(GDALGetRasterBand(m_dataset, i + 1), GF_Read, 0, 0, xsize, ysize, imgBuffer, xsize, ysize, GDT_UInt16, 0, 0);
		for (int j = 0; j < ysize; j++)
		{
			for (int k = 0; k < xsize; k++)
			{
				//扩大100倍精度为0.01
				imgBuffer[k] = (unsigned short)((imgBuffer[j*xsize + k] * parametersA[i] + parametersB[i]) * 100);
			}
		}
		GDALRasterIO(GDALGetRasterBand(m_datasetdst, i + 1), GF_Write, 0, 0, xsize, ysize, imgBuffer, xsize, ysize, GDT_UInt16, 0, 0);
	}

	delete[]parametersA; parametersA = NULL;
	delete[]parametersB; parametersB = NULL;
	delete[]imgBuffer;	 imgBuffer = NULL;
	GDALClose(m_dataset);
	GDALClose(m_datasetdst);
	return lError;
}

GDALDataset* geGdalVSI::VsiGdalCreateWrap(GDALDriverH hdriver,
                                          std::string* const vsifile,
                                          int nx, int ny, int nbands,
                                          GDALDataType bandtype,
                                          char **papszoptions) {
  GDALDatasetH hvsi_ds;
  khMutex &mutex = GetMutex();
  khLockGuard lock(mutex);
  *vsifile = UniqueVSIFilename();
  hvsi_ds = GDALCreate(hdriver, (*vsifile).c_str(), nx, ny, nbands,
                       bandtype, papszoptions);
  return static_cast<GDALDataset*>(hvsi_ds);
}

GDALDatasetH QgsRelief::openOutputFile( GDALDatasetH inputDataset, GDALDriverH outputDriver )
{
  if ( !inputDataset )
  {
    return nullptr;
  }

  int xSize = GDALGetRasterXSize( inputDataset );
  int ySize = GDALGetRasterYSize( inputDataset );

  //open output file
  char **papszOptions = nullptr;

  //use PACKBITS compression for tiffs by default
  papszOptions = CSLSetNameValue( papszOptions, "COMPRESS", "PACKBITS" );

  //create three band raster (reg, green, blue)
  GDALDatasetH outputDataset = GDALCreate( outputDriver, mOutputFile.toUtf8().constData(), xSize, ySize, 3, GDT_Byte, papszOptions );
  if ( !outputDataset )
  {
    return outputDataset;
  }

  //get geotransform from inputDataset
  double geotransform[6];
  if ( GDALGetGeoTransform( inputDataset, geotransform ) != CE_None )
  {
    GDALClose( outputDataset );
    return nullptr;
  }
  GDALSetGeoTransform( outputDataset, geotransform );

  //make sure mCellSizeX and mCellSizeY are always > 0
  mCellSizeX = geotransform[1];
  if ( mCellSizeX < 0 )
  {
    mCellSizeX = -mCellSizeX;
  }
  mCellSizeY = geotransform[5];
  if ( mCellSizeY < 0 )
  {
    mCellSizeY = -mCellSizeY;
  }

  const char *projection = GDALGetProjectionRef( inputDataset );
  GDALSetProjection( outputDataset, projection );

  return outputDataset;
}

GDALDatasetH QgsRasterCalculator::openOutputFile( GDALDriverH outputDriver )
{
  //open output file
  char **papszOptions = NULL;
  GDALDatasetH outputDataset = GDALCreate( outputDriver, mOutputFile.toLocal8Bit().data(), mNumOutputColumns, mNumOutputRows, 1, GDT_Float32, papszOptions );
  if ( outputDataset == NULL )
  {
    return outputDataset;
  }

  //assign georef information
  double geotransform[6];
  outputGeoTransform( geotransform );
  GDALSetGeoTransform( outputDataset, geotransform );

  return outputDataset;
}

GDALDatasetH QgsRasterCalculator::openOutputFile( GDALDriverH outputDriver )
{
  //open output file
  char **papszOptions = nullptr;
  GDALDatasetH outputDataset = GDALCreate( outputDriver, TO8F( mOutputFile ), mNumOutputColumns, mNumOutputRows, 1, GDT_Float32, papszOptions );
  if ( !outputDataset )
  {
    return outputDataset;
  }

  //assign georef information
  double geotransform[6];
  outputGeoTransform( geotransform );
  GDALSetGeoTransform( outputDataset, geotransform );

  return outputDataset;
}

gdal::dataset_unique_ptr QgsRasterCalculator::openOutputFile( GDALDriverH outputDriver )
{
  //open output file
  char **papszOptions = nullptr;
  gdal::dataset_unique_ptr outputDataset( GDALCreate( outputDriver, mOutputFile.toUtf8().constData(), mNumOutputColumns, mNumOutputRows, 1, GDT_Float32, papszOptions ) );
  if ( !outputDataset )
  {
    return nullptr;
  }

  //assign georef information
  double geotransform[6];
  outputGeoTransform( geotransform );
  GDALSetGeoTransform( outputDataset.get(), geotransform );

  return outputDataset;
}

/* Makes a copy of a dataset, and opens it for writing.. */
GDALDatasetH make_me_a_sandwitch(GDALDatasetH *in_dataset, char *filename)
{
    char **papszOptions = NULL;
    const char *pszFormat = "GTiff";
    GDALDriverH hDriver;
    GDALDatasetH out_gdalfile;
    hDriver = GDALGetDriverByName( pszFormat );
    papszOptions = CSLSetNameValue( papszOptions, "TILED", "YES" );
    papszOptions = CSLSetNameValue( papszOptions, "COMPRESS", "DEFLATE" );
    
    /*Create copy..*/
    /*return GDALCreateCopy( hDriver, filename, *in_dataset, FALSE, papszOptions, NULL, NULL );*/
    return GDALCreate(hDriver, filename,
        GDALGetRasterXSize( *in_dataset ),
        GDALGetRasterYSize( *in_dataset ),
        1,
        GDT_Byte, papszOptions );
}

static dErr JakoGDALDatasetCreateMem(OGRSpatialReferenceH ref,const double geo[6],dInt n,dInt nlines,GDALDataType dtype,GDALDatasetH *dset,void *bandmem)
{
  char *wkt;
  GDALDriverH memdriver;
  CPLErr cplerr;
  OGRErr oerr;
  dErr err;

  dFunctionBegin;
  oerr = OSRExportToWkt(ref,&wkt);dOGRCHK(oerr);
  memdriver = GDALGetDriverByName("MEM");
  *dset = GDALCreate(memdriver,"MEM:::",n,nlines,0,dtype,NULL);
  cplerr = GDALSetProjection(*dset,wkt);dCPLCHK(cplerr);
  cplerr = GDALSetGeoTransform(*dset,(double*)geo);dCPLCHK(cplerr); /* const-incorrect interface */
  OGRFree(wkt);
  if (bandmem) {err = JakoGDALMemAddBand(*dset,GDT_Float64,&bandmem);dCHK(err);}
  dFunctionReturn(0);
}

s_sat * sat_rasterize_copy(s_sat * c_sat, OGRGeometryH geometry)
{
	try;

		int band = 1;
		unsigned height, width, height_width;
		double burn_value = 255;
		GDALDatasetH ds = NULL;
		s_sat * sat = NULL;

		throw_null((sat = sat_init(1)));

		height = sat->height = c_sat->height;
		width = sat->width = c_sat->width;
		height_width = height * width;
		sat->proj_ref = strdup(c_sat->proj_ref);
		sat->ch_num = 0;

		memcpy(sat->gt_coef, c_sat->gt_coef, sizeof(double) * 6);
		throw_null((sat->pixel[0] = sfire_alloc(sizeof(uint8_t), 1, height_width)));

		sat->ch_num = 1;

		throw_null((ds = GDALCreate(drv_r_mem, "", width, height, 1, GDT_Byte, NULL)));
		throw((GDALSetProjection(ds, sat->proj_ref) == CE_Failure));
		throw((GDALSetGeoTransform(ds, sat->gt_coef) == CE_Failure));

		throw((GDALRasterizeGeometries(ds, 1, & band, 1, & geometry, NULL, NULL, & burn_value, NULL, NULL, NULL) == CE_Failure));

		throw(GDALDatasetRasterIO(ds, GF_Read, 0, 0, width, height, sat->pixel[0], width, height, GDT_Byte, 1, NULL, 0, 0, 0) == CE_Failure);

	catch;

		sat_destroy(sat);

		sat = NULL;

	finally;

		if(ds != NULL)
			GDALClose(ds);

	return sat;
}

GDALDatasetH QgsNineCellFilter::openOutputFile( GDALDatasetH inputDataset, GDALDriverH outputDriver )
{
  if ( inputDataset == NULL )
  {
    return NULL;
  }

  int xSize = GDALGetRasterXSize( inputDataset );
  int ySize = GDALGetRasterYSize( inputDataset );;

  //open output file
  char **papszOptions = NULL;
  GDALDatasetH outputDataset = GDALCreate( outputDriver, mOutputFile.toLocal8Bit().data(), xSize, ySize, 1, GDT_Float32, papszOptions );
  if ( outputDataset == NULL )
  {
    return outputDataset;
  }

  //get geotransform from inputDataset
  double geotransform[6];
  if ( GDALGetGeoTransform( inputDataset, geotransform ) != CE_None )
  {
    GDALClose( outputDataset );
    return NULL;
  }
  GDALSetGeoTransform( outputDataset, geotransform );

  //make sure mCellSizeX and mCellSizeY are always > 0
  mCellSizeX = geotransform[1];
  if ( mCellSizeX < 0 )
  {
    mCellSizeX = -mCellSizeX;
  }
  mCellSizeY = geotransform[5];
  if ( mCellSizeY < 0 )
  {
    mCellSizeY = -mCellSizeY;
  }

  const char* projection = GDALGetProjectionRef( inputDataset );
  GDALSetProjection( outputDataset, projection );

  return outputDataset;
}