c_cpp 遥感图像拼块处理框架

Posted 2021-05-14

#include <gdal_priv.h>
#include <Eigen\Dense>

int numY = 0;
int lineCount = height;
if ( width > 1000 ) { lineCount = 50; }
if ( width > 10000 ) { lineCount = 5 ; }
int currentLineCount = lineCount;
int tileCount = height / lineCount;

for ( numY = 0; numY < tileCount + 1; numY++ )
{
    // deal with last piece of tile
    if ( numY == tileCount )
    {
        currentLineCount = height % lineCount;
        if ( currentLineCount == 0 ) { break; }
    }

    // read data out
    float** srcData = new float*[this->bandCount];
    MatrixXd src( bandCount, this->width * currentLineCount );
    for ( int band = 0; band < this->bandCount; band++ )
    {
        srcData[band] = new float[this->width * currentLineCount];
        int bandList = {band + 1};
        inputDataset->RasterIO( GF_Read, 0, numY * lineCount, this->width, currentLineCount, srcData[band],
                                width, currentLineCount, GDT_Float32, 1, &bandList, 0, 0, 0 );
        int noDataValue;
        GDALGetRasterNoDataValue( inputDataset->GetRasterBand( bandList ), &noDataValue );
        // x
        for ( int j = 0; j < this->width * currentLineCount; j++ )
        {
            if( srcData[band][j] != noDataValue &&
                    _finite( srcData[band][j] ) ) { src( band, j ) = srcData[band][j]; }
            else { src( band, j ) = 0; }
        }
    }
    delete []srcData;

    float* res = new float[width * currentLineCount]; // result array
    
    // TO DO: write processing code here
    res = ;
    
    
    pand->RasterIO( GF_Write, 0, numY * lineCount, this->width, currentLineCount, res,
                    this->width, currentLineCount, GDT_Float32, 0, 0 );
}

void writeMatOut( std::string matName, Eigen::MatrixXd &mat )
{
    std::string dir = "";
    ofstream out( dir + matName + ".txt", ios::out );
    if ( !out.is_open() ) { throw std::string( "write file error." ); }

    out << mat ;
    out.close();
}