Project

General

Profile

o2-osmocom-bb-smscb2kml.pl

roox, 02/07/2018 08:29 PM

 
1
#!/usr/bin/perl
2

    
3
#############################################################################
4
## >>> 2014 <<<<
5
## Martin Hauke (mardnh@gmx.de)
6
##
7
## Older GSM BTS's from O2 Germany (formerly MNC 07 - now MNC 03) send out their BTS
8
## coordinates via smscb (channel 221).
9
## With Alex Badea's patches for libosmocore/osmocombb applied you can log those smscb
10
## messages via osmocom-bb's cell_log application.
11

    
12
## This script parses osmocom-bb cell_log output and writes the relevant info to a
13
## nice kml-file.
14
##
15
## usage example:
16
## cat cell.log | o2-osmocom-bb-smscb2kml.pl > nice_map.kml
17
##
18
##
19
## The used Converter for Gauss - Krueger coordinates is licensed under the GPL and was
20
## written by
21
## - Andreas Achtzehn (achtzehn@linux-society.de)
22
## - Nobert Hüttisch (nobbi@nobbi.com)
23

    
24

    
25
#use strict;
26

    
27
my %seen;
28

    
29
my ($line, $time, $arfcn, $mcc, $mnc, $lac, $cid, $cb);
30
my ($viag_gk, $sr, $sx);
31
my $regex = qw(^CellBroadcast: TIME=(?<time>\d+)\sARFCN=(?<arfcn>\d+)\sMCC=(?<mcc>\d+)\sMNC=(?<mnc>\d+)\sLAC=(?<lac>.*)\sCID=(?<cid>.*)\sCB=\"(?<cb>.*)\"$);
32

    
33
&printKMLheader;
34

    
35
foreach $line (<STDIN>) {
36
	chomp($line);
37
	if ( $line =~ /$regex/ ) {
38
		#print "$line\n";
39
		$time  = "$+{time}";
40
		$arfcn = "$+{arfcn}";
41
		$mcc   = "$+{mcc}";
42
		$mnc   = "$+{mnc}";
43
		$lac   = "$+{lac}";
44
		$cid   = "$+{cid}";
45
		$cb    = "$+{cb}";
46
		$cid_dec = hex($cid);
47
		$lac_dec = hex($lac);
48
		$time_x  = scalar localtime($time);
49

    
50
		# check if O2-Germany (O2 MNC was until mid 2016 - since then 03)
51
		if ($mnc == "03" && $mcc == "262") {
52
			# did we seen that entry point already?
53
			next if $seen{"$arfcn,$mcc,$mnc,$lac,$cid"}++;
54
			$cb =~ /(?<gk>\d+).*/;
55
			$viag_gk = "$+{gk}";
56
			$xlat=substr($viag_gk,0,6)*10;
57
			$xlon=substr($viag_gk,6,6)*10;
58
			getCoordinate($xlat, $xlon);
59
			#print "Google Maps: https://maps.google.com/maps?z=12&h=m&q=loc:$wgs84lat+$wgs84lon&z=15\n";
60
			#print "$time $arfcn $mcc $mnc $lac $cid $viag_gk $wgs84lat $wgs84lon\n";
61
			
62
			# print placemarks
63
			print qq{\t\t<Placemark>\n};
64
			print qq{\t\t\t<name>CID:$cid_dec</name>\n};
65
			print qq{\t\t\t<description><![CDATA[CID: <b>$cid ($cid_dec)</b><br> LAC: <b>$lac ($lac_dec)</b><br>ARFCN: <b>$arfcn</b><br> MCC: <b>$mcc</b><br>MNC: <b>$mnc</b><br>TIME: <b>$time_x</b> ]]></description>\n};
66
			#print qq{\t\t\t</description>\n};
67
			#print "<styleUrl>#rxl69</styleUrl>"
68
			print qq{\t\t\t<Point>\n};
69
			print qq{\t\t\t\t<coordinates>$wgs84lon,$wgs84lat</coordinates>\n};
70
			print qq{\t\t\t</Point>\n};
71
			print qq{\t\t</Placemark>\n};
72
		};
73
	};
74
};
75

    
76
&printKMLfooter;
77

    
78
sub getCoordinate {
79
	$sr = shift;
80
	$sx = shift;
81
	&potsdam;
82
	&modolensky;
83
}
84

    
85
sub potsdam {
86
	# Errechnen der Position nach Potsdam-Koordinatensystem
87
	$bm  = int($sr/1000000);
88
	$y   = $sr-($bm*1000000+500000);
89
	$si  = $sx/111120.6196;
90
	$px  = $si+0.143885358*sin(2*$si*0.017453292)+0.00021079*sin(4*$si*0.017453292)+0.000000423*sin(6*$si*0.017453292);
91
	$t   = (sin($px*0.017453292))/(cos($px*0.017453292));
92
	$v   = sqrt(1+0.006719219*cos($px*0.017453292)*cos($px*0.017453292));
93
	$ys  = ($y*$v)/6398786.85;
94
	$lat = $px-$ys*$ys*57.29577*$t*$v*$v*(0.5-$ys*$ys*(4.97-3*$t*$t)/24);
95
	$dl  = $ys*57.29577/cos($px*0.017453292) * (1-$ys*$ys/6*($v*$v+2*$t*$t-$ys*$ys*(0.6+1.1*$t*$t)*(0.6+1.1*$t*$t)));
96
	$lon = $bm*3+$dl;
97

    
98
	$ph =int($lat);
99
	$pz =($lat-$ph)*60;
100
	$pm =int($pz);
101
	$ps =($pz-$pm)*60;
102

    
103
	$lh =int($lon);
104
	$lz =($lon-$lh)*60;
105
	$lm =int($lz);
106
	$ls =($lz-$lm)*60;
107
}
108

    
109
sub modolensky {
110
	$potsd_a  = 6377397.155;
111
	$wgs84_a  = 6378137.0;
112
	$potsd_f  = 1/299.152812838;
113
	$wgs84_f  = 1/298.257223563;
114

    
115
	$potsd_es = 2*$potsd_f - $potsd_f*$potsd_f;
116

    
117
	$potsd_dx = 606.0;
118
	$potsd_dy = 23.0;
119
	$potsd_dz = 413.0;
120
	$pi = 3.141592654;
121
	$latr = $lat/180*$pi;
122
	$lonr = $lon/180*$pi;
123

    
124
	$sa = sin($latr);
125
	$ca = cos($latr);
126
	$so = sin($lonr);
127
	$co = cos($lonr);
128

    
129
	$bda  = 1-$potsd_f;
130

    
131
	$delta_a = $wgs84_a - $potsd_a;
132
	$delta_f = $wgs84_f - $potsd_f;
133

    
134
	$rn = $potsd_a / sqrt(1-$potsd_es*sin($latr)*sin($latr));
135
	$rm = $potsd_a * ((1-$potsd_es)/sqrt(1-$potsd_es*sin($latr)*sin($latr)*1-$potsd_es*sin($latr)*sin($latr)*1-$potsd_es*sin($latr)*sin($latr)));
136

    
137
	$ta = (-$potsd_dx*$sa*$co - $potsd_dy*$sa*$so)+$potsd_dz*$ca;
138
	$tb = $delta_a*(($rn*$potsd_es*$sa*$ca)/$potsd_a);
139
	$tc = $delta_f*($rm/$bda+$rn*$bda)*$sa*$ca;
140
	$dlat = ($ta+$tb+$tc)/$rm;
141

    
142
	$dlon = (-$potsd_dx*$so + $potsd_dy*$co)/($rn*$ca);
143

    
144
	$wgs84lat = ($latr + $dlat)*180/$pi;
145
	$wgs84lon = ($lonr + $dlon)*180/$pi;
146
	$ph=int($wgs84lat);
147
	$pz=($wgs84lat-$ph)*60;
148
	$pm=int($pz);
149
	$ps=int(($pz-$pm)*60);
150
	
151
	$lh=int($wgs84lon);
152
	$lz=($wgs84lon-$lh)*60;
153
	$lm=int($lz);
154
	$ls=int(($lz-$lm)*60);
155

    
156
	return ($wgs84lat, $wgs84lon )
157
}
158

    
159
sub printKMLheader {
160
	print qq{<?xml version="1.0" encoding="UTF-8"?>\n};
161
	print qq{<kml xmlns="http://www.opengis.net/kml/2.2">\n};
162
	print qq{\t<Document>\n};
163
};
164

    
165
sub printKMLfooter {
166
	print qq{\t</Document>\n};
167
	print qq{</kml>\n};
168
};
169

    
170
# EOF
Add picture from clipboard (Maximum size: 48.8 MB)