source: subversion/applications/utils/osm-matrix/osm-matrix.pl @ 34613

Last change on this file since 34613 was 22068, checked in by tappenbeck, 9 years ago

osm-matrix public

File size: 8.0 KB
Line 
1#####################################################################################################
2# OSM-MATRIX - create a gpx-grid for show the grid of osm-matrix by monty                                                                                                                                #
3#                                                                                                                                                                                                                                                       #
4# Copyright (C) 2010 Jan Tappenbeck, osm(at)tappenbeck.net                                                                                                                                                              #
5# components in use of garry68 (gpx-file-componetents), components for tile-calculation (http://wiki.openstreetmap.org/wiki/Slippy_map_tilenames)                               #
6#                                                                                                                                                                                                                                                       #
7# This program is free software; you can redistribute it and/or modify it under the terms of the                                                                                                        #
8# GNU General Public License as published by the Free Software Foundation; either version 3 of                                                                                                  #
9# the License, or (at your option) any later version.                                                                                                                                                                           #
10# This program is distributed in the hope that it will be useful, but WITHOUT ANY WARRANTY;                                                                                                     #
11# without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.                                                                                     #
12# See the GNU General Public License for more details.                                                                                                                                                                  #
13# You should have received a copy of the GNU General Public License along with this program;                                                                                                            #
14# if not, see <http://www.gnu.org/licenses/>.                                                                                                                                                                                   #       
15#                                                                                                                                                                                                                                                       #
16#  DEVELOP and TESTED in WINDOWS VISTA 64bit / ActivePerl                                                                                                                                               #
17#                                                                                                                                                                                                                               #
18#####################################################################################################
19
20#*****************
21# HISTORIE
22#*****************
23# 2010-01-25  first in use by haiti-earthquake
24
25#!/usr/bin/perl
26# permalink für test
27
28use strict;
29use LWP;
30#use OSM::osm;
31#use OSM::osm_jt ; #paketbasis von gary68
32use Math::Trig;
33use Getopt::Long;
34
35my $go_Help                     = 0;            # help please
36my $matrix_zoom         = 15;
37my $counter             = 0;
38my $output_prefix       = "";
39my $coord_west          = -9999;
40my $coord_north         = -9999;
41my $coord_east          = -9999;
42my $coord_south         = -9999;
43
44GetOptions
45(
46        "help!"                                 => \$go_Help,
47        "name=s"                                => \$output_prefix,
48        "w=f"                                   => \$coord_west,
49        "n=f"                                   => \$coord_north,
50        "e=f"                                   => \$coord_east,
51        "s=f"                                   => \$coord_south,
52        "z=i"                                   => \$matrix_zoom,
53) or Usage();
54
55Usage() if( $go_Help );                                                                 # Hilfe wenn erwuenscht
56
57#----------------------------------------------------------------------------------------
58sub Usage
59{
60        my( $message ) = @_;
61
62        if( $message )
63        {
64                print "$message\n";
65        }
66
67        print "\tgpx-grid for osm-matrix\n\n";
68       
69        print "\tparameters\n";
70       
71        print "\tosm [s]\tinputfile\n";
72        print "\tname [s]\tprefix for the gpx-file -> osm_matrix-[name].gpx\n";
73        print "\n\tmin/max geographic coordinates\n";
74        print "\tw [f]\twest-limit\n";
75        print "\tn [f]\tnorth-limit\n";
76        print "\te [f]\teast-limit\n";
77        print "\ts [f]\tsouth-limit\n";
78
79        print "\n\tz [i]\ttile-zoom - default 15\n";
80       
81        print "\n\n\t-- end --\n";
82
83        exit 0;                                 # beenden
84}
85
86
87my $output_gpx_filename = "osm_matrix_".$output_prefix.".gpx";
88
89if ($coord_east < $coord_west){
90  print "*"x51,"\n";
91  print "******* the area change the date-border !!! *******\n";
92  print "******* please split the area into 2 grid *******\n";
93  print "*"x51,"\n";
94  exit 0;
95}
96
97#if ($coord_north < $coord_south){
98#  print "*"x59,"\n";
99#  print "******* north < south - please check the parameters *******\n";
100#  print "*"x59,"\n";
101#  exit 0;
102#}
103
104# transform geo. coord into tile-coord. for the define zoom-scale
105
106
107my ($lat_start,$lon_end) = getTileNumber($coord_south, $coord_west, $matrix_zoom);
108my ($lat_end,$lon_start) = getTileNumber($coord_north, $coord_east, $matrix_zoom);
109
110# to get the right border of the selected area
111$lon_end++;
112
113print "\n\ngeo-coord.\n";
114print "buttom-left: lon=".$coord_south."  lat=".$coord_west."\n";
115print "top-right  : lon=".$coord_north." lat=".$coord_east."\n\n";
116
117print "\n\ntile-numbers\n";
118print "buttom-left: lat=".$lat_start." / lon=".$lon_start."\n";
119print "top-right  : lat=".$lat_end." / lon=".$lon_end."\n\n";
120
121open (my $gpxfile, ">", $output_gpx_filename) || die ("Can't open html output file: ".$output_gpx_filename) ;
122printGPXHeader($gpxfile);
123
124print "==> horizontal lines\n";
125 for( my $h=$lon_start; $h <= $lon_end; $h++ )
126   {
127        print $lat_start." / ".$h." - ".$lat_end." / ".$h."\n";
128        my ($lat1, $lon1, $lat_tmp1, $lon_tmp1) =  Project($lat_start,$h,$matrix_zoom);
129        my ($lat2, $lon2, $lat_tmp2, $lon_tmp2) =  Project($lat_end,$h,$matrix_zoom);
130    print $lat1." / ".$lon1." - ".$lat2." / ".$lon2."\n";
131        printOpenGpsTrack($gpxfile, $counter);
132        printGPXTrackpoint($gpxfile, $lon1, $lat1);
133        printGPXTrackpoint($gpxfile, $lon2, $lat2);
134        printCloseGpsTrack($gpxfile);
135  }
136
137print "==> vertical lines\n";
138for( my $r=$lat_start; $r <= $lat_end; $r++ )
139{
140    print $r." / ".$lon_start." - ".$r." / ".$lon_end."\n";
141        my ($lat1, $lon1, $lat_tmp1, $lon_tmp1) =  Project($r,$lon_start,$matrix_zoom);
142        my ($lat2, $lon2, $lat_tmp2, $lon_tmp2) =  Project($r,$lon_end,$matrix_zoom);
143    print $lat1." / ".$lon1." - ".$lat2." / ".$lon2."\n";
144        printOpenGpsTrack($gpxfile, $counter);
145        printGPXTrackpoint($gpxfile, $lon1, $lat1);
146        printGPXTrackpoint($gpxfile, $lon2, $lat2);
147        printCloseGpsTrack($gpxfile);
148} 
149
150printGPXFoot($gpxfile);
151
152#*****************
153# SUBROUTINES
154#*****************
155
156############
157# tile-calculation
158############
159
160# source for the calculation of the request coord.
161# http://projekte.eiops.de/osm-matrix/?zoom=16&lat=18.49757&lon=-72.63565&layers=B0FTFFFFFFFFFT
162#  http://wiki.openstreetmap.org/wiki/Slippy_map_tilenames
163
164sub getTileNumber {
165  my ($lat,$lon,$zoom) = @_;
166  my $xtile = int( ($lon+180)/360 *2**$zoom ) ;
167  my $ytile = int( (1 - log(tan($lat*pi/180) + sec($lat*pi/180))/pi)/2 *2**$zoom ) ;
168  return(($xtile, $ytile));
169}
170
171 sub Project {
172  my ($X,$Y, $Zoom) = @_;
173  my $Unit = 1 / (2 ** $Zoom);
174  my $relY1 = $Y * $Unit;
175  my $relY2 = $relY1 + $Unit;
176
177  # note: $LimitY = ProjectF(degrees(atan(sinh(pi)))) = log(sinh(pi)+cosh(pi)) = pi
178  # note: degrees(atan(sinh(pi))) = 85.051128..
179  #my $LimitY = ProjectF(85.0511);
180
181  # so stay simple and more accurate
182  my $LimitY = pi;
183  my $RangeY = 2 * $LimitY;
184  $relY1 = $LimitY - $RangeY * $relY1;
185  $relY2 = $LimitY - $RangeY * $relY2;
186  my $Lat1 = ProjectMercToLat($relY1);
187  my $Lat2 = ProjectMercToLat($relY2);
188  $Unit = 360 / (2 ** $Zoom);
189  my $Long1 = -180 + $X * $Unit;
190  return(($Lat2, $Long1, $Lat1, $Long1 + $Unit)); # S,W,N,E
191 }
192 sub ProjectMercToLat($){
193  my $MercY = shift();
194  return( 180/pi* atan(sinh($MercY)));
195 }
196 sub ProjectF
197 {
198  my $Lat = shift;
199  $Lat = deg2rad($Lat);
200  my $Y = log(tan($Lat) + (1/cos($Lat)));
201  return($Y);
202 }
203
204########
205# gpx-track
206########
207# legt einen GPS-Track an
208sub printGPXHeader {
209        my $file = shift ;
210
211        print $file "<?xml version=\"1.0\" encoding=\"UTF-8\" standalone=\"no\" ?>\n" ;
212        print $file "<gpx xmlns=\"http://www.topografix.com/GPX/1/1\" creator=\"Gary68script\" version=\"1.1\"\n" ;
213        print $file "    xmlns:xsi=\"http://www.w3.org/2001/XMLSchema-instance\"\n" ;
214        print $file "    xsi:schemaLocation=\"http://www.topografix.com/GPX/1/1 http://www.topografix.com/GPX/1/1/gpx.xsd\">\n" ;
215}
216sub printGPXFoot {
217        my $file = shift ;
218
219        print $file "</gpx>\n" ;
220}
221sub printGPXWaypoint {
222        my ($file, $lon, $lat, $text) = @_ ;
223
224        print $file "<wpt lat=\"", $lat, "\" lon=\"", $lon, "\">" ;
225        print $file "<desc>", $text, "</desc></wpt>\n" ;
226}
227sub printOpenGpsTrack {
228        my $file = shift ;
229        my $name = shift ;
230
231        print $file " <trk>\n";
232        print $file "  <name>$name</name>\n";
233        print $file "   <trkseg>\n";
234}
235sub printCloseGpsTrack {
236        my $file = shift ;
237
238        print $file "  </trkseg>\n" ;
239        print $file "</trk>\n" ;
240}
241sub printGPXTrackpoint {
242        my ($file, $lon, $lat) = @_ ;
243        print $file "    <trkpt lat=\"", $lat, "\" lon=\"", $lon, "\">\n";
244        print $file "    </trkpt>\n" ;
245}
Note: See TracBrowser for help on using the repository browser.