Using gdal_translate for multiband raster to convert into .XYZ file?

Question

I have raster containing 2 bands (NDVI and DEM). This raster was created in SNAP using "Add Elevation Band" option, "DEM Model:SRTM 1sec Grid" to the NDVI product I had generated earlier. I could not find any option in SNAP to make the conversion of this combined product into an ASCII file. There is however one option in SNAP to export the product to a CSV but the result of that process were not as expected. So I switched to QGIS. Here I tried the Raster>Conversion>Translate option, which ultimately uses the gdal_translate tool.

So in QGIS-Translate I could convert the .tif to .XYZ(ascii) but the output had only one band value along with the lat-long for each pixel. I got the warning "Warning 6: XYZ driver only uses the first band of the dataset."

I used the option -b 1 -b 3 to include the bands I want in my output file but still the same warning kept coming.

My gdal_query looks as below,

gdal_translate -of XYZ test-op-1p-ndvi_DEM.tif gdaltranslate-op1DEM.xyz -b 1 -b 3

I tried using the other option in this query like -ot with Byte format (as mentioned in one related question here in GIS-StackExchange) and also by removing the -of XYZ parameter

I have also tried gdal_translate through the OSGeo4W desktop tool but got same warning there too,

After searching hard enough, I could not find anything satisfactory. So now I am trying the long way by generating 2 separate .XYZ for each band and then combine through python script, parsing and combining the 2 columns.

Can anyone please provide any hint or suggestion. Even if someone can link to similar questions in this forum would be really helpful, as i have already tried solutions from many related questions but could have missed something important.

Screenshot of .XYZ with only NDVI values

Screenshot of .XYZ with only DEM values

Answer 1

GDAL XYZ driver supports only single band data. The documentation https://www.gdal.org/frmt_xyz.html somehow suggests that more bands (more than one value per cell) could be supported

Those datasets are ASCII files with (at least) 3 columns, each line containing the X and Y coordinates of the center of the cell and the value of the cell.

However, by readind the source code https://github.com/OSGeo/gdal/blob/master/gdal/frmts/xyz/xyzdataset.cpp it appears that when reading the data other columns than those recognized as X, Y, and Z are skipped. When writing exactly one source band is required:

/* -------------------------------------------------------------------- */
/*      Some some rudimentary checks                                    */
/* -------------------------------------------------------------------- */
    int nBands = poSrcDS->GetRasterCount();
    if (nBands == 0)
    {
        CPLError( CE_Failure, CPLE_NotSupported,
                  "XYZ driver does not support source dataset with zero band.\n");
        return nullptr;
    }

    if (nBands != 1)
    {
        CPLError( (bStrict) ? CE_Failure : CE_Warning, CPLE_NotSupported,
                  "XYZ driver only uses the first band of the dataset.\n");
        if (bStrict)
            return nullptr;
    }

    if( pfnProgress && !pfnProgress( 0.0, nullptr, pProgressData ) )
return nullptr;

However, there is a python utility https://github.com/OSGeo/gdal/blob/master/gdal/swig/python/scripts/gdal2xyz.py that you can use instead of gdal_translate and zyx driver. Read more from Exporting DEM-derived rasters to CSV from QGIS?.