wiki:PerlMapScriptExamples35ex1

create_tileindex.pl

#!/usr/bin/perl
use strict;
#use warnings;

use XBase;
use Cwd;
use mapscript;
use Getopt::Long;
use Image::Size;
use File::Basename;
use File::Find;
use Pod::Usage;

our $version = 0.91;
our $ms37above = 0;             # default mapserver version < 3.7

BEGIN {
  my ($msmajor, $msminor) = ($mapscript::MS_VERSION =~ /^(\d+)\.(\d+)/o);
  if ($msmajor > 3 or ($msmajor == 3 and $msminor > 6)) {
    $ms37above = 1;
  }
}

# variables for command line arguments
my $verbose   = 0;              # verbosity level
my $imagedir  = getcwd;         # image directory, default is current dir
my $shapename = 'tileindex';    # name of shapefile to be created
my $extension = 'tif';          # extension of image files
my $imgcatdrv = '';             # image catalog drive and path
my $imgcatsbd = '';             # image catalog share basepath
my $help      = '';             # you want help? you'll get it!

GetOptions('verbose+'        => \$verbose,
           'imagedir:s'      => \$imagedir,
           'tileindexname:s' => \$shapename,
           'extension:s'     => \$extension,
           'drvimgcat:s'     => \$imgcatdrv,
           'sbdimgcat:s'     => \$imgcatsbd,
           'help'            => \$help);

pod2usage( { -exitval => 1,
             -verbose => 3} ) if $help;

pod2usage( { -exitval => 2,
             -verbose => 1,
             -output => \*STDERR} ) unless $imagedir and $shapename and $extensi
on;

my $NL = "\n";

# If tileindex shape exists, delete it!
for my $shpfile ("$shapename.shp", "$shapename.shx", "$shapename.dbf") {
  if (-e $shpfile) {
    print "deleting $shpfile\n" if $verbose > 1;
    unlink $shpfile
      or die "ERROR: Could not delete $shpfile: $!\n";
  }
}

die "ERROR: Not a directory $imagedir: $!\n" unless -d $imagedir;
my @tfws;
find( \&push_tfws, $imagedir);
die "ERROR: No .tfw files in directory $imagedir\n" if $#tfws < 0;

# create new shapefile
print "creating tileindex $shapename.shp\n\n" if $verbose > 1;
my $shapefile = $ms37above ? mapscript::shapefileObj->new("$shapename.shp",$mapscript::MS_SHAPEFILE_POLYGON) : shapefileObj->new("$shapename.shp",$mapscript::MS_SHAPEFILE_POLYGON);
my $basename = $shapename =~ /\.shp$/ ? basename($shapename,('.shp')) : $shapename;
# create dBASE-file for shapefile
unlink "$shapename.dbf";
my $table =  XBase->create("name"        => "$shapename.dbf",
                  "field_names"    => [ qw/ IMAGE XMIN YMIN XMAX YMAX LOCATION /
 ],
                  "field_types"    => [ qw/ C N N N N C/ ],
                  "field_lengths"  => [ qw/ 75 20 20 20 20 127/ ],
                  "field_decimals" => [ undef, 8, 8, 8, 8, undef ]);

# create point object for usage in shapefile
my $point = $ms37above ? mapscript::pointObj->new() : pointObj->new() ;

my $i = 0;                      # record counter
foreach my $tfwfile (@tfws) {   # step through all tfw-files
  my $image = imagename($tfwfile); # name of image file
  print "Image is $image...\n" if $verbose;
  next unless -r $image;        # skip if image does not exist
  print "Processing $tfwfile...\n" if $verbose;
  # calculate image coordinates
  my $tfw = read_tfw("$tfwfile", $image);
  # create new shape object which will contain the new polygon
  my $shp = $ms37above ? mapscript::shapeObj->new($mapscript::MS_POLYGON) : shapeObj->new($mapscript::MS_POLYGON);
  # create new line object which will receive the coordinates
  my $line = $ms37above  ? mapscript::lineObj->new() : lineObj->new();
  my ($minx, $miny, $maxx, $maxy) = ();
  foreach my $rk (@$tfw) {      # add all coordinates
    print $$rk[0], ',', $$rk[1], $NL if $verbose > 1;
    $point->{'x'} = $$rk[0];    # into point object
    $point->{'y'} = $$rk[1];
    $line->add($point);         # and point into line
    save_maxmin($rk, \$minx, \$miny, \$maxx, \$maxy);
  }
  # write record into dbase file
  (my $avimage = $image) =~ s!$imgcatsbd!$imgcatdrv!o;
  $avimage =~ s!/!\\!og;
  $table->set_record($i,$avimage,$minx,$miny,$maxx,$maxy,$image);
  $shp->add($line);             # add line to shape
  $shapefile->add($shp);        # add shape to shapefile
  $i++;                         # increase record counter
  print $NL if $verbose > 1;
}

undef $shapefile;               # close shapefile
undef $table;                   # close dbase table
system('shptree', $shapename);
print "Done.\n" if $verbose;
exit 0;

#---------------------------------------------------------------------
# function save_maxmin - save maximum/minimum coordinates
#
# description: compare current point koordinates with max/min x/y values
#              and set these new if appropriate
#
# parameters : $rk   - reference to coordinate list
#              $minx - reference to min x value
#              $miny - reference to min y value
#              $maxx - reference to max x value
#              $maxy - reference to max y value
#
# returns    : nothing
#
# MWS remarks:
#
#---------------------------------------------------------------------
sub save_maxmin {
  my ($rk, $minx, $miny, $maxx, $maxy) = @_;
  $$minx = $$minx ? $$rk[0] < $$minx ? $$rk[0] : $$minx : $$rk[0];
  $$maxx = $$maxx ? $$rk[0] > $$maxx ? $$rk[0] : $$maxx : $$rk[0];
  $$miny = $$miny ? $$rk[1] < $$miny ? $$rk[1] : $$miny : $$rk[1];
  $$maxy = $$maxy ? $$rk[1] > $$maxy ? $$rk[1] : $$maxy : $$rk[1];
  print STDERR "$$minx, $$miny, $$maxx, $$maxy\n" if $verbose > 1;
}


#---------------------------------------------------------------------
# function imagename - return name of image for given tfw file
#
# description: removes extension .tfw and adds known extension
#              from global variable $extension
#
# parameters : $img - name of tfw file
#
# returns    : string with image name
#---------------------------------------------------------------------
sub imagename{
  my $img = shift;
  $img =~ s/\.tfw$//o;
  return "$img.$extension";
}

#---------------------------------------------------------------------
# function read_tfw - read tfwfile and image, return coordinates
#
# description: reads tfw file and retrieves image size, creates a list
#              of the corner coordinates
#              Uses Image::Size to get the size of the image
#
# parameters : $tfwfile - name of tfw file
#              $image   - name of image
#
# returns    : pointer to list of coordinates
#
#---------------------------------------------------------------------
sub read_tfw{
  my ($tfwfile, $image) = @_;
  my ($ix, $iy, $desc) = imgsize($image);
  die "$desc" unless $ix;
  open(I, "< $tfwfile") or die "Can't read $tfwfile: $!";
  my @lines = <I>;              # read file into array
  my ($dx) =  ($lines[0] =~ /(\S+)/)[0];
  my ($dy) =  ($lines[3] =~ /(\S+)/)[0];
  my ($ulx) = ($lines[4] =~ /(\S+)/)[0];
  my ($uly) = ($lines[5] =~ /(\S+)/)[0];
  close I;
  my $coords = [[$ulx,$uly],
                  [$ulx + $dx * $ix, $uly],
                  [$ulx + $dx * $ix, $uly + $dy * $iy],
                  [$ulx, $uly + $dy * $iy],
                  [$ulx,$uly]];
  return $coords;
}

#---------------------------------------------------------------------
# function push_tfws - push all tfw-files on global stack @tfw
#
# description:
#
# parameters :
#
# returns    :
#
# MWS remarks:
#
#---------------------------------------------------------------------
sub push_tfws{
  push(@tfws, $File::Find::name) if /\.tfw$/;
}

1;

__END__

NAME

create_tileindex.pl - create tileindex for mapserver and image catalog for ArcView

SYNOPSIS

  create_tileindex.pl [--tileindexname hitif] [--extension tif] [--imagedir /tmp/hitif] [--drvimgcat y:\]                        
                      [--sbdimgcat /data] [--verbose] [--help]

DESCRIPTION

create_tileindex.pl creates a tileindex shapefile and dbase file for the usage of tiled georeferenced images in mapserver and ArcView? (TM ESRI).

The images can be stored in a subdirectory tree, as F<create_tileindex.pl> performs a file-find starting from the given imagedir.

After creating the tileindex, a F<shptree> command is issued to speed up mapserver access to the tiles.

EXAMPLE

Simple usage

Cd into base directory of the images, issue the F<create_tileindex.pl> command. A tileindex named F<tileindex> will be created in the current directory. This is sufficient for mapserver usage and ArcView? usage with a F<drivemap.txt> file.

Complex usage

  create_tileindex.pl --tileindexname muctiles --extension tif --imagedir /data/tifplan --sbdimgcat /data                          
                          --drvimgcat "y:"

Short notation:

create_tileindex.pl -t muctiles -i /data/tifplan -s /data -d "y:"

This command creates the tileindex with the name F<muctiles> in the current directory.

The extension of the searched image files is tif (this is the default value, so it could have been omitted).

The starting directory of the search for image files is F</data/tifplan>.

For the usage as an image catalog, the directory part of the image filenames is modified by

=over 4

=item

replacing from the start the value of the --sbdimgcat value with the --drvimgcat value

=item

replacing all forward slashes / with backslashes =back

The modified filenames are placed in the item C<IMAGE>. It's length is 75 characters.

(For ArcView? experts: The same effect could be achieved by placing an appropriate F<DRIVEMAP.TXT> in the appropriate place.)

=head1 COMMAND LINE ARGUMENTS

Every argument name can be given in one-dash-one-letter notation or shortened ad libitum.

=over 4

=item --tileindexname (optional)

Name (without extension!) of the tileindex to be created. If the tileindexname contains no directory part, the tileindex is created in the current directory. If this argument is omitted, the name F<tileindex> is used.

=item --extension (optional)

Extension of the image files to be searched. Defaults to C<tif>.

=item --imagedir (optional)

Base directory where the image files are stored under. The search for the images starts here.

=item drvimgcat (optional)

Drive letter and eventually path (Windows notation) to be used when accessing the images from ArcView/Windows? via the image catalog.

=item sbdimgcat (optional)

Windows share base directory of the image files (Unix notation). This directory will be replaced by the C<drvimgcat> value.

=item verbose (optional)

Toggles verbose output. One -v argument prints the names of the processed image files, two or more additionally print the corner coordinates of the images.

=item help (optional)

Print this documentation to stdout and exit.

=back

Tileindex in MapServer

F<create_tileindex.pl> places the names of the image file in the item C<LOCATION>, which is the default for mapserver. It's length is 127 characters. To use the created tileindex in mapserver, add the following lines in your mapfile:

  LAYER
    NAME bgl99
    TILEINDEX "/home/springm/perl/muctiles"
    TYPE raster
  END

Image catalog in ArcView

The created dbase file can be used as an image catalog in ArcView?. To use this image catalog, perform the following steps (quoted from ArcView? help):

=over 4

=item

Click the Add Theme button.

=item

In the Data Source Types box, choose Image Data Source.

=item

Navigate to the directory that contains the image catalog you want to add. Double-click on the directory name to list the image catalogs sources it contiains.

=item

Double-click the image catalog source you wish to add.

=back

=head1 LIMITATIONS

F<create_tileindex.pl> does not yet work with geotifs a.k.a. tif files which include the reference information in their header.

=head1 REQUIREMENTS

F<create_tileindex.pl> uses perl 5.6.1 and is based on C<mapscript>, Version 3.6 or above.

It further uses the Perl modules C<XBase>, C<Image::Size> and <Pod::Usage>. You can get them from L<CPAN|http://search.cpan.org>.

=head1 AUTHOR

Markus W. Spring <m.spring@…>

=head1 COPYRIGHT AND DISCLAIMER

This program is Copyright 2002.2003 by Markus W. Spring. This program is free software; you can redistribute it and/or modify it under the terms of the GNU General Public License as published by the Free Software Foundation; either version 2 of the License, or (at your option) any later version.

This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more details.

If you do not have a copy of the GNU General Public License write to the Free Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.

=cut

# Local variables: # compile-command: "perl create_tileindex.pl -t ihk_tk50 -i /data/ihk/tk50 -e tif -v -v" # End:


back to PerlMapScript

Last modified 13 years ago Last modified on Jan 28, 2009, 2:04:14 PM