pyrsgis 0.3.2

A package to read, process and export GeoTIFFs.

Python-for-Remote-Sensing-and-GIS

pyrsgis enables the user to read, process and export GeoTIFFs. The module is built on the GDAL library, but is much more convenient when it comes to reading and exporting GeoTIFs. pyrsgis also supports reading satellite data directly from .tar.gz files. However, reading from .tar.gz files is currently in its beta phase. Please do not use this package for commercial purpose without my explicit permission. Researchers/ academicians are welcomed for feedback and technical support. Since this is an open-source volunatry project, collaborations are most welcome. Please write to me at [e-mail]

See installation command using pip on the PyPi page - link

Recommended citation:
Tripathy, P. pyrsgis: A Python package for remote sensing and GIS. V0.3.2

Sample code

1. Reading .tif extension file

Import the module and define the input file path.

from pyrsgis import raster

file_path = r'D:/your_file_name.tif'

• To read all the bands of a stacked satellite image:
  

ds, arr = raster.read(file_path, bands='all')

where, ds is the data source similar to GDAL and arr is the numpy array that contains all the bands of the input raster. The arr can be 2D or 3D depending on the input data. One can check the shape of the array using the print(arr.shape) command. The bands argument in the raster.read function defaults to 'all'.

• To read a list of bands of a stacked satellite image:
  

ds, arr = raster.read(file_path, bands=[2, 3, 4])

Passing the band numbers in a list returns bands 2, 3 & 4 as three-dimensional numpy array.

• To read a specific band from stacked satellite image:
  

ds, arr = raster.read(file_path, bands=2)

Passing a single band number returns that particular band as two-dimensional numpy array.

• To read a single band TIF file:
  

ds, arr = raster.read(file_path)

Since the bands argument defaults to 'all', this will read all the bands in the input file, here, one band only.

2. Exporting .tif extension file

In all the below examples, it is assumed that the number of rows and columns, and the cell size of the input and output rasters are the same. All these are stored in the `ds` variable, please see details here: link.

• To export all bands of a 3D array:
  

out_file_path = r'D:/sample_file_all_bands.tif'
raster.export(arr, ds, out_file_path, dtype='int', bands='all')

The dtype argument in the above function defaults to 'int', which is 8-bit integer. Please be careful to change this while exporting arrays with large values. Similarly, to export float type array (eg. NDVI), use dtype = 'float'. Data type of high pixel-depth, e.g. Integer32, Integer64, or float type uses more space on hard drive, so the default has been set to integer. To export any float datatype, the argument should be passed explicitly.
These are the options for the dtype argument: 'byte', 'cfloat32', 'cfloat64', 'cint16', 'cint32', 'float32', 'float64', 'int16', 'int32', 'uint8', 'uint16', 'uint32'

• To export a list of bands of a 3D array:
  

out_file_path = r'D:/sample_file_bands_234.tif'
raster.export(arr, ds, out_file_path, bands=[2, 3, 4])

• To export any one band of a 3D array:
  

out_file_path = r'D:/sample_file_band_3.tif'
raster.export(arr, ds, out_file_path, bands=3)

• To export a single band array:
  

out_file_path = r'D:/sample_file.tif'
raster.export(arr, ds, out_file_path)

where, arr should be a 2D array.

3. Converting TIF to CSV

GeoTIFF files can be converted to CSV files using *pyrsgis*. Every band is flattened to a single-dimensional array, and converted to CSV. These are very useful for statistical analysis.
Import the function:

from pyrsgis.convert import rastertocsv

To convert all the bands present in a folder:

your_dir = r"D:/your_raster_directory"
out_file_path = r"D:/yourFilename.csv"

rastertocsv(your_dir, filename=out_file_path)

Generally the NoData or NULL values in the raster become random negative values, negatives can be removed using the negative argument:

rastertocsv(your_dir, filename=out_file_path, negative=False)

At times the NoData or NULL values in raster become '127' or '65536', they can also be removed by declaring explicitly.

rastertocsv(your_dir, filename=out_file_path, remove=[127, 65536])

This is a trial and check process, please check the generated CSV file for such issues and handle as required.

Similarly, there are bad rows in the CSV file, representing zero value in all the bands. This takes a lot of unnecessary space on drive, it can be eliminated using:

rastertocsv(your_dir, filename=out_file_path, badrows=False)

4. Creating northing and easting using a reference raster

pyrsgis allows to quickly create the northing and easting rasters using a reference raster, as shown below:

The cell value in the generated rasters are the row and column number of the cell. To generate these GeoTIFF files, start by importing the function:

from pyrsgis.raster import northing, easting

reference_file_path = r'D:/your_reference_raster.tif'

northing(reference_file_path, outFile= r'D:/pyrsgis_northing.tif', flip=True)
easting(reference_file_path, outFile= r'D:/pyrsgis_easting.tif', flip=False)

As the name suggests, the flip argument flips the resulting rasters.

5. Reading directly from .tar.gz files (beta)

Currently, only Landsat data is supported.

import pyrsgis

file_path = r'D:/your_file_name.tar.gz'
your_data = pyrsgis.readtar(file_path)

The above code reads the data and stores in the your_data variable.

Various properties of the raster can be assessed using the following code:

print(your_data.rows)
print(your_data.cols)

This will display the number of rows and columns of the input data.

Similarly, the number of bands can be checked using:

print(your_data.nbands)

On reading the .tar.gz files directly, pyrsgis determines the satellite sensor. This can be checked using:

print(your_data.satellite)

This will display the satellite sensor, for instance, Landsat-5, Landsat-8, etc.

If the above code shows the correct satellite sensor, then the list of band names of the sensor (in order) can easily be checked using:

print(your_data.bandIndex)

Any particular band can be extarcted using:

band_number = 1
your_band = your_data.getband(band_number)

The above code returns the band as array which can be visualised using:

display(your_band, maptitle='Title of your image', cmap='PRGn')

or, directly using:

band_number = 1
display(your_data.getband(band_number), maptitle='Title of your image', cmap='PRGn')

The generated map can directly be saved as an image.

The extracted band can be exported using:

out_file_path = r'D:/sample_output.tif'
your_data.export(your_band, out_file_path)

This saves the extracted band to the same directory.

To export the float type raster, please define the datatype explicitly, default is 'int':

your_data.export(your_band, out_file_path, datatype='float')

The NDVI (Normalised Difference Vegetaton Index) can be computed easily.

your_ndvi = your_data.ndvi()

Normalised difference index between any two bands can be computed using:

norm_diff = your_data.nordif(bandNumber2, bandNumber1)

This computes (band2-band1)/(band2+band1) in the back end and returns a numpy array. THe resulting arracy can be exported using:

out_file_path = r'D:/your_ndvi.tif'
your_data.export(your_ndvi, out_file_path, datatype='float')

Be careful with the float type of NDVI.