3

I am using MATLAB in my office and Octave when I am at home. Although they are very similar, I was trying to do something I would expected to be very easy and obvious, but found it really annoying. I can't find out how to import TIFF images in Octave. I know the MATLAB geotiffread function is not present, but I thought there would be another method.

I could also skip importing them, as I can work with the imread function in some cases, but then the second problem would be that I can't find a way to write a georeferenced TIFF file (in MATLAB I normally call geotiffwrite with geotiffinfo inputs inside). My TIFF files are usually 8 bit unsigned integer or 32 bit signed integer. I hope someone can suggest a way to solve this problem. I also saw this thread but did not understand if it is possible to use the code proposed by Ashish in Octave.

Cris Luengo
  • 55,762
  • 10
  • 62
  • 120
umbe1987
  • 2,894
  • 6
  • 35
  • 63

3 Answers3

3

You may want to look at the mapping library in Octave. You can also use the raster functions to work with GeoTiffs

Example:

pkg load mapping
filename=”C:\\sl\SDK\\DTED\\n45_w122_1arc_v2.tif”
rasterinfo (filename)
rasterdraw (filename)
2

The short answer is you can't do it in Octave out of the box. But this is not because it is impossible to do it. It is simply because no one has yet bothered to implement it. As a piece of free software, Octave has the features that its users are willing to spend time or money implementing.

About writing of signed 32-bit images

As of version 3.8.1, Octave uses either GraphicsMagick or ImageMagick to handle the reading and writing of images. This introduces some problems. The number 1 is that your precision is limited to how you built GraphicsMagick (its quantum-depth option). In addition, you can only write unsigned integers. Hopefully this will change in the future but since not many users require it, it's been this way until now.

Dealing with geotiff

Provided you know C++, you can write this functions yourself. This shouldn't be too hard since there is already libgeotiff, a C library for it. You would only need to write a wrapper as an Octave oct function (of course, if you don't know C or C++, then this "only" becomes a lot of work).

carandraug
  • 12,938
  • 1
  • 26
  • 38
0

Here is the example oct file code which needs to be compiled. I have taken reference of https://gerasimosmichalitsianos.wordpress.com/2018/01/08/178/

#include <octave/oct.h>
#include "iostream"
#include "fstream"
#include "string"
#include "cstdlib"
#include <cstdio>

#include "gdal_priv.h"
#include "cpl_conv.h"
#include "limits.h"
#include "stdlib.h"
using namespace std;
typedef std::string String;

DEFUN_DLD (test1, args, , "write geotiff")
{
  NDArray maindata = args(0).array_value ();
  const dim_vector dims = maindata.dims ();
  int i,j,nrows,ncols;
    nrows=dims(0);
    ncols=dims(1);  
        //octave_stdout << maindata(i,0);
    NDArray transform1 = args(1).array_value ();
  double*  transform = (double*) CPLMalloc(sizeof(double)*6);
  float*  rowBuff = (float*) CPLMalloc(sizeof(float)*ncols);
  //GDT_Float32 *rowBuff = CPLMalloc(sizeof(GDT_Float32)*ncols);
  String tiffname;
  tiffname  = "nameoftiff2.tif";
  
  cout<<"The transformation matrix is";
  for (i=0; i<6; i++)
  {
      transform[i]=transform1(i);
      cout<<transform[i]<<" ";
  }
     
  GDALAllRegister();
  CPLPushErrorHandler(CPLQuietErrorHandler);
  GDALDataset *geotiffDataset;
  GDALDriver *driverGeotiff;
  GDALRasterBand *geotiffBand;
  
  OGRSpatialReference oSRS;
  char **papszOptions = NULL;
  char *pszWKT = NULL;
 
    oSRS.SetWellKnownGeogCS( "WGS84" );
    oSRS.exportToWkt( &pszWKT );
    
  driverGeotiff = GetGDALDriverManager()->GetDriverByName("GTiff");  
  geotiffDataset = (GDALDataset *) driverGeotiff->Create(tiffname.c_str(),ncols,nrows,1,GDT_Float32,NULL);
  geotiffDataset->SetGeoTransform(transform);
  geotiffDataset->SetProjection(pszWKT);
  //CPLFree( pszSRS_WKT );
  cout<<" \n Number of rows and columns in array are: \n";
  cout<<nrows<<" "<<ncols<<"\n";
  for (i=0; i<nrows; i++)
    {
        for (j=0; j <ncols; j++)
            rowBuff[j]=maindata(i,j);
        //cout<<rowBuff[0]<<"\n";
        geotiffDataset->GetRasterBand(1)->RasterIO(GF_Write,0,i,ncols,1,rowBuff,ncols,1,GDT_Float32,0,0);
    }
  
  
  GDALClose(geotiffDataset) ;
  CPLFree(transform);
  CPLFree(rowBuff);
  CPLFree(pszWKT);
  GDALDestroyDriverManager();
  
  return octave_value_list();
}

it can be compiled and run using following

mkoctfile -lgdal test1.cc
aa=rand(50,53);
b=[60,1,0,40,0,-1];
test1(aa,b);