#!/usr/bin/python
import sys
from osgeo import gdal
try:
    raster = gdal.Open(sys.argv[1])
    gt = raster.GetGeoTransform()
    k = abs(gt[1] * gt[5] - gt[4] * gt[2])
    band = raster.GetRasterBand(int(sys.argv[2]))
    nodata = band.GetNoDataValue() 
    r = band.ReadAsArray()
    if nodata is not None:
        r = r[r != nodata]
    if len(sys.argv) > 3:
        r = r[r >= float(sys.argv[3])]
    if len(sys.argv) > 4:
        r = r[r < float(sys.argv[4])]
    print(k * r.size)
except IndexError:
    print('Usage: {0} RASTER_FILE BAND_NO [MIN [MAX]]'.format(sys.argv[0]))
    print('Returns total area of raster pixels in range [MIN, MAX)')