package Geo::GSHHS;
use strict;
use Carp;


sub new()
  {
    my $self  = {};
    $self->{headerformat}   = 
      "N" . # Polygon ID
	"N" . # Number of points
	  "N" . # Level (1 = land)
	    "N4" . # West, east, north, south
	      "N" . # Area in tenths of kilometers-squared
		"n" . # 1 if Greenwich meridian is crossed
		  "n";  # Source data
    $self->{pointdataformat} = "NN"; # X, Y
    $self->{polygons} = ();
    $self->{numpoints} = 0;
    $self->{numpolygons} = 0;
    $self->{next} = 0;

    bless($self);
    return $self;

}

sub get_numpolygons()
  {
    my $self = shift();
    return($self->{numpolygons});
  }

sub get_numpoints()
  {
    my $self = shift();
    return($self->{numpoints});
  }

sub next
  {
    my $self = shift();
    return(0) if ($self->{next} >= $self->{numpolygons});
    my $polygon = ${$self->{polygons}}[$self->{next}];
    $self->{next}++;

    return($polygon);
  }

sub reset()
  {
    my $self = shift();
    $self->{next} = 0;
    return(1);
  }

sub close()
  {
    my $self = shift();
    $self->{polygons} = ();
    return(1);
  }

#----------------------------------------------
# Opens and reads a GSHHS file
# returns true if it succeeds
#----------------------------------------------
sub open()
  {
    my $self = shift();
    my $filename = shift();
    my @polygons = ();
    my $totalpoints = 0;
    my $totalpolygons = 0;

    # Open the file for reading, binary
    open(my $fp, "<:bytes", $filename) || return(0);

    # Continue reading headers until end of file
    while(read($fp, my $header, 8 * 4 + 2 * 2))
      {

	# Read the binary data into local variables
	my ($id, 
	    $numpoints,
	    $level,
	    $west,
	    $east,
	    $north,
	    $south,
	    $area,
	    $crosses,
	    $source) = unpack($self->{'headerformat'}, $header);


	# New blank array, local
	my @points = ();

	# Loop through vertices in this polygon
	foreach(1..$numpoints)
	  {
	    
	    # Read from file
	    if(read($fp, my $datapoint, 2 * 4))
	      {

		# Unpack binary data into local variables
		my($x, $y) = unpack($self->{'pointdataformat'}, $datapoint);

		my $latitude = sign($y) / 1E+6;
		my $longitude = sign($x) / 1E+6;
		$longitude -= 360 if($longitude > 180);

		my $point = {longitude => $longitude,
			     latitude => $latitude};

		# Add lat/long as a hash to the array of vertices
		push(@points, $point);
		
	      }
	    else
	      {
		carp "Error in GSHHS data: too many vertices in a polygon\n";
	      }

	  }

	# Create a local hash using the header variables, and point data
	my $polygon = {'numpoints' => $numpoints,
		       level     => sign($level),
		       west      => sign($west)   / 1E+6,
		       east      => sign($east)   / 1E+6,
		       north     => sign($north)  / 1E+6,
		       south     => sign($south)  / 1E+6,
		       area      => $area   / 10,
		       crosses   => $crosses,
		       points    => \@points};


	# Push onto the list of polygons
	push(@polygons, $polygon);
	$totalpoints += $numpoints;
	$totalpolygons++;

      }

    # Close the file
    close($fp);

    # Store the list of polygons
    $self->{polygons} = \@polygons;
    $self->{numpolygons} = $totalpolygons;
    $self->{numpoints} = $totalpoints;

    # Return success
    return(1);
  }


1;

sub sign()
  {
    my $x = shift();
    if($x > 0x80000000)
      {
	$x = $x - 0xFFFFFFFF;
	$x--;
      }
    return($x);
  }

__END__

=head1 NAME

Geo::GSHHS

A module to read Global Self-consistant Hierarchical High-resolution 
Shorelines data, see ftp://gmt.soest.hawaii.edu/pub/pwessel/gshhs/
for details, and for the data files.


=head1 Functions


=over 4

=item new()

Creates a new instance of the parser

=item open(filename)

Reads a GSHHS binary file from disk, and stores
it internally.

 Arguments  : filename

 Returns: 1 : succeeded
          0 : failed

=item close()

Deletes any internal data being stored, which may be quite large.

=item get_numpolygons()

Returns the number of polygons in the internal data

=item get_numpoints()

Returns the number of vertices in the internal data

=item reset()

Resets an internal  pointer, so that next() will start
at the beginning of the data set

=item next()

Retrieves the next polygon from the data. (see data description,
below)


=back

=head1 Data


 - numpolygons
 - polygons
   -   numpoints
   -   level:         1 land,
	              2 lake,
	              3 island_in_lake, 
	              4 pond_in_island_in_lake
   -   west           westernmost limit in degrees
   -   east           easternmost
   -   north          nothernmost
   -   south          southernmost
   -   area           area in kilometers-squared
   -   crosses        spans greenwich meridian?
   -   points[numpoints]
        -  latitude   degrees (north positive)
        -  longitude  degrees (east positive)


=head1 Example
 use Geo::GSHHS;

 my $Coastlines = Geo::GSHHS::new();
 $Coastlines->open('data/gshhs_i.b');
 $Coastlines->reset();
 while(my $Polygon = $Coastlines->next())
 {
   foreach my $count(0 .. $Polygon->{numpoints} - 1)
     {
       print $Polygon->{points}->[$count]->{latitude}, "\t",
	 $Polygon->{points}->[$count]->{longitude}, "\n";

     }
 }

=head1 AUTHOR

Oliver White <oliver.white@blibbleblobble.co.uk>.


=cut back
