File: mapblast2geo.pl

package info (click to toggle)
xastir 2.1.4-1
  • links: PTS, VCS
  • area: main
  • in suites: bullseye, sid
  • size: 13,676 kB
  • sloc: ansic: 121,100; perl: 7,529; sh: 5,535; makefile: 383; python: 312; sql: 102
file content (199 lines) | stat: -rwxr-xr-x 7,270 bytes parent folder | download | duplicates (2)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
#!/usr/bin/perl

# XASTIR .geo file generator for mapblast pixel maps    16.10.2001
#  Copyright (C) 2001 Rolf Bleher              http://www.dk7in.de

#  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.
#
#  You should have received a copy of the GNU General Public License
#  along with this program; if not, write to the Free Software
#  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA

#  see file COPYING for details

#--------------------------------------------------------------------------

use POSIX;                      # provides acos() function
                                            
# some default values:
$width  = 1280;                 # seems to be the maximum supported by mapblast
$height = 1024;                 # seems to be the maximum supported by mapblast
$pfx    = "../../pixelmaps/";   # prefix for path
$pfx    = "";

# undefined values:
$scale  =    0;
$lat    = 9999;
$lon    = 9999;
$help   =    0;

#--------------------------------------------------------------------------

while (@ARGV) {
    $arg = shift @ARGV;
    if ($arg =~ /^\-N(\d{1,2}(\.\d+)?)$/)  { $lat   =  $1 }
    if ($arg =~ /^\-S(\d{1,2}(\.\d+)?)$/)  { $lat   = -$1 }
    if ($arg =~ /^\-E(\d{1,3}(\.\d+)?)$/)  { $lon   =  $1 }
    if ($arg =~ /^\-W(\d{1,3}(\.\d+)?)$/)  { $lon   = -$1 }
    if ($arg =~ /^\-w(\d+)$/)              { $width =  $1 }
    if ($arg =~ /^\-h(\d+)$/)              { $height=  $1 }
    if ($arg =~ /^\-s(\d+)$/)              { $scale =  $1 }
    if ($arg =~ /^\-s(\d+)k$/)             { $scale =  $1 * 1000 }
    if ($arg =~ /^\-s(\d+)M$/)             { $scale =  $1 * 1000000 }
    if ($arg =~ /^\-p(.+)$/)               { $pfx   =  $1 }
    if ($arg eq '-?')                      { $help  =  1 }
    if ($arg eq '-h')                      { $help  =  1 }
    if ($arg eq '--help')                  { $help  =  1 }
}

if (!$help && $arg && $arg !~ /^\-/) {
    $file = $arg;                       # last arg is file name
} else {
    print("ERROR: Map file name is missing\n");
    $help = 1;
}

if (!$help) {
    if ($lat == 9999) {         # we don't yet have a latitude
        if ($file =~ /N(\d{1,2}(\.\d+)?)/)  { $lat =  $1 }
        if ($file =~ /S(\d{1,2}(\.\d+)?)/)  { $lat = -$1 }
    }
    if ($lon == 9999) {         # we don't yet have a longitude
        if ($file =~ /E(\d{1,3}(\.\d+)?)/)  { $lon =  $1 }
        if ($file =~ /W(\d{1,3}(\.\d+)?)/)  { $lon = -$1 }
    }
    if ($scale == 0) {          # we don't yet have a map scale
        if ($file =~ /\-(\d+)/)  { $scale =  $1 }
        if ($file =~ /\-(\d+)k/) { $scale =  $1 * 1000 }
        if ($file =~ /\-(\d+)M/) { $scale =  $1 * 1000000 }
    }
    if ($lat == 9999) {         # we need a latitude
        print("ERROR: Latitude is missing\n");
        $help = 1;
    }
    if ($lon == 9999) {         # we need a longitude
        print("ERROR: Longitude is missing\n");
        $help = 1;
    }
    if ($scale == 0) {          # we need a map scale
        print("ERROR: Map scale is missing\n");
        $help = 1;
    }
}
    
if ($help) {
    usage();
    exit 0;
}

if ($pfx && $pfx !~ /\/$/) { $pfx .= "/" }      # add trailing '/' if missing

#--------------------------------------------------------------------------

if ($scale >= 1000000) {
    printf("ERROR: Maps with scaling of 1M and above are not yet supported!\n");
    printf("       They use a different projection...\n");
    exit 1;
}

if (abs($lat) > 89) {
    printf("ERROR: Map center too near to the poles!\n");
    exit 1;
}

# This scaling factor is just a wild guess!
# But I need one, and this one is nice AND works quite well... ;-)
$scale_y = 100.0 * pi();                # pixel/degree at 1M map scale

$scale_x = $scale_y * calcScale($lat);  # adjust horizontal scale for latitude

$scale_y *= (1000000 / $scale);         # adjust for current map scale
$scale_x *= (1000000 / $scale);

#--------------------------------------------------------------------------

# DK7IN: I'm not sure, if this formula is exact for what Xastir
#        is decoding, but in my region it works quite well
#        I need to do some further investigation for the best accuracy
$latmin = $lat - ($height / 2.0 - 2.5) / $scale_y;
$latmax = $lat + ($height / 2.0 - 0.5) / $scale_y;
$lonmin = $lon - ($width  / 2.0 - 0.5) / $scale_x;
$lonmax = $lon + ($width  / 2.0 - 2.5) / $scale_x;

printf("FileName %s%s\n",$pfx,$file);
printf("\n");
printf("ImageSize  %4d %4d\n",$width,$height);
printf("TiePoint   %4d %4d  %10.6f %11.6f\n",0,0,$lonmin,$latmax);
printf("TiePoint   %4d %4d  %10.6f %11.6f\n",$width-1,$height-1,$lonmax,$latmin);
printf("Datum      WGS84\n");
printf("Projection LatLon\n");
printf("\n");
printf("# map from mapblast center %10.6f %11.6f, scale %d\n",$lat,$lon,$scale);
printf("# created with mapblast2geo.pl (DK7IN)\n");

exit 0;

#--------------------------------------------------------------------------

sub usage {
    my $name = $0;
    if ($name =~ /^.*\/(.+)$/) { $name = $1 }
    print("\n");
    print("$name (c) 2001 Rolf Bleher <Rolf\@dk7in.de>\n");
    print("create Xastir .geo files for mapblast pixel maps\n");
    print("usage: $name [options] mapfile\n");
    print("       -N52.5 -S10    define latitude\n");
    print("       -E13.3 -W0.5   define longitude\n");
    print("       -h1024         define map height in pixels (default 1024)\n");
    print("       -w1280         define map width in pixels  (default 1280)\n");
    print("       -s50000        define map scale, 50k or 1M is ok\n");
    print("       -p../pixmaps   define prefix for path\n");
    print("       -h -? --help   print this help file\n");
    print("    it tries to extract center Lat/Lon and map scale\n");
    print("    from the filename like N52.5E13.3-50k.xpm\n");
    print("\n");
}

#--------------------------------------------------------------------------

sub pi {
    return(3.14159265358979323846);
}

#--------------------------------------------------------------------------

sub deg2rad {
    my ($deg) = @_;
    return($deg * pi()/180.0);
}

#--------------------------------------------------------------------------

# Calculate distance in meters between two locations
sub dist {
    my ($lat1, $lon1, $lat2, $lon2) = @_;
    my $r_lat1 = deg2rad($lat1);
    my $r_lon1 = deg2rad($lon1);
    my $r_lat2 = deg2rad($lat2);
    my $r_lon2 = deg2rad($lon2);
    my $r_d = acos(sin($r_lat1) * sin($r_lat2) + cos($r_lat1) * cos($r_lat2) * cos($r_lon1-$r_lon2));
    return($r_d*180.0*60.0/pi()*1852.0);
}

#--------------------------------------------------------------------------

sub calcScale {                 # EW / NS scaling
    my ($lat) = @_;
    return(dist($lat,-1.0/120.0,$lat,1.0/120.0) / 1852.0);
}

#--------------------------------------------------------------------------