Geo::GDAL
1.11
|
The code here does not have any warranty. It is recommended that before using any of this code, you look into it and try to understand what it does, what input it needs, etc. Do not blindly execute anything!
use Carp; use Geo::GDAL; use PDL; my $filename = 'xxx'; my $update = undef; my $dataset = Geo::GDAL::Open($filename, $update); my $band = $dataset->Band(1); my($xoff,$yoff) = (1000,4000); my($w,$h) = (1000,1000); my $data = $band->ReadRaster($xoff,$yoff,$w,$h); my $pdl = PDL->new; # PDL data types: # 0 = unsigned char # 1 = short # 2 = unsigned short # 3 = int # 4 = long long # 5 = float # 6 = double # GDAL data types: Byte UInt16 Int16 UInt32 Int32 Float32 Float64 CInt16 CInt32 CFloat32 CFloat64 # obtained with $band->DataType my %map = ( Byte => 0, UInt16 => 2, Int16 => 1, UInt32 => -1, Int32 => 3, Float32 => 5, Float64 => 6, CInt16 => -1, CInt32 => -1, CFloat32 => -1, CFloat64 => -1 ); my $datatype = $map{$band->DataType}; croak "there is no direct mapping between the band datatype and PDL" if $datatype < 0; $pdl->set_datatype($datatype); $pdl->setdims([1,$w,$h]); my $dref = $pdl->get_dataref(); # this is not optimal, but the best what we can do in Perl: $$dref = $data; $pdl->upd_data;