This repository has been archived by the owner on May 7, 2018. It is now read-only.
forked from oscarpilote/Ortho4XP
-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathOrtho4XP.py
8552 lines (8140 loc) · 379 KB
/
Ortho4XP.py
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
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
#!/usr/bin/env python3
##############################################################################
# Ortho4XP : A base mesh creation tool for the X-Plane 11 flight simulator. #
# Version : 1.20b #
# Copyright 2016 Oscar Pilote #
# Thanks to all that have contributed to improvement of the code. #
##############################################################################
# #
# LEGAL NOTICE : #
# #
# 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 3 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, see <http://www.gnu.org/licenses/>. #
# #
##############################################################################
version=' v1.20b'
import os
import sys
try:
import encodings.idna
except:
pass
if getattr(sys,'frozen',False):
Ortho4XP_dir = '..'
os.environ["REQUESTS_CA_BUNDLE"] = os.path.join(os.getcwd(), "cacert.pem") # needed to access https providers with the bin version
else:
Ortho4XP_dir = '.'
try:
sys.path.append(os.getcwd()+'/'+Ortho4XP_dir+'/bin/Modules')
except:
pass
import requests
import threading,subprocess,time,gc,shutil,io
from math import *
import array,numpy
import random
import collections
import struct
import hashlib
from tkinter import * # GUI
from tkinter import filedialog
import tkinter.ttk as ttk # Themed Widgets
from PIL import Image, ImageDraw, ImageFilter, ImageTk
Image.MAX_IMAGE_PIXELS = 1000000000 # Not a decompression bomb attack!
import subprocess
try:
import gdal
gdal_loaded = True
except:
gdal_loaded = False
########################################################################
#
# FACTORY DEFAULT VALUES
#
# The following are initialisation of the variables, they are then superseded by your Ortho4XP.cfg file but are here
# in case your config file misses some of then (evolution, edit, ...). Do not modify them here but use your
# config file instead.
configvars=['default_website','default_zl','water_option','sea_texture_params','cover_airports_with_highres',\
'cover_zl','cover_extent','min_area','sea_equiv','do_not_flatten_these_list','meshzl',\
'insert_custom_zoom_in_mesh','poly_simplification_tol','overpass_server_list','overpass_server_choice',\
#
'curvature_tol','no_small_angles','smallest_angle','hmin','hmax','tile_has_water_airport',\
'water_smoothing','sea_smoothing',\
#
'masks_width','complex_masks','use_masks_for_inland','legacy_masks','keep_old_pre_mask',\
'maskszl','use_gimp','gimp_cmd',\
#
'ratio_water','normal_map_strength','terrain_casts_shadows','use_decal_on_terrain',\
'dds_or_png','use_bing_for_non_existent_data', 'max_convert_slots','be_nice_timer',\
'skip_downloads','skip_converts','check_tms_response',\
'verbose_output','clean_tmp_files','clean_unused_dds_and_ter_files',\
'contrast_adjust','brightness_adjust','saturation_adjust',\
'full_color_correction','g2xpl_8_prefix','g2xpl_8_suffix','g2xpl_16_prefix','g2xpl_16_suffix',\
#
'Custom_scenery_prefix','Custom_scenery_dir','default_sniff_dir',\
'keep_orig_zuv','seven_zip','landclass_mesh_division','snap_to_z_grid','overlay_lod',\
'keep_overlays'\
]
configvars_strings=['default_website','overpass_server_choice','gimp_cmd','dds_or_png',\
'g2xpl_8_prefix','g2xpl_8_suffix','g2xpl_16_prefix','g2xpl_16_suffix',\
'Custom_scenery_prefix','Custom_scenery_dir','default_sniff_dir']
configvars_defaults={\
'default_website':'BI',\
'default_zl':16,\
'water_option':3,\
'sea_texture_params':[],\
'cover_airports_with_highres':False,\
'cover_zl':18,\
'cover_extent':1,\
'min_area':0.01,\
'sea_equiv':[],\
'do_not_flatten_these_list':[],\
'meshzl':19,\
'insert_custom_zoom_in_mesh':False,\
'poly_simplification_tol':0.002,\
'overpass_server_list':{"1":"http://api.openstreetmap.fr/oapi/interpreter", "2":"http://overpass-api.de/api/interpreter","3":"http://overpass.osm.rambler.ru/cgi/interpreter"},\
'overpass_server_choice':1,\
'curvature_tol':3,\
'no_small_angles':False,\
'smallest_angle':5,\
'hmin':20,\
'hmax':2000,\
'tile_has_water_airport':False,\
'water_smoothing':2,\
'sea_smoothing':0,\
'masks_width':16,\
'complex_masks':False,\
'use_masks_for_inland':False,\
'legacy_masks':True,\
'keep_old_pre_mask':False,\
'maskszl':14,\
'use_gimp':False,\
'gimp_cmd':'',\
'ratio_water':0.3,\
'normal_map_strength':0.3,\
'terrain_casts_shadows':False,\
'use_decal_on_terrain':False,\
'dds_or_png':'dds',\
'use_bing_for_non_existent_data':False,\
'max_convert_slots':4,\
'be_nice_timer':0,\
'skip_downloads':False,\
'skip_converts':False,\
'check_tms_response':True,\
'verbose_output':True,\
'clean_tmp_files':True,\
'clean_unused_dds_and_ter_files':False,\
'contrast_adjust':{},\
'brightness_adjust':{},\
'saturation_adjust':{},\
'full_color_correction':{},\
'g2xpl_8_prefix':'g2xpl_8_',\
'g2xpl_8_suffix':'',\
'g2xpl_16_prefix':'g2xpl_16_',\
'g2xpl_16_suffix':'',\
'Custom_scenery_prefix':'',\
'Custom_scenery_dir':'',\
'default_sniff_dir':'',\
'keep_orig_zuv':True,\
'seven_zip':True,\
'landclass_mesh_division':8,\
'snap_to_z_grid':True,\
'overlay_lod':40000,\
'keep_overlays':True\
}
explanation={}
for item in configvars:
explanation[item]='TODO!'
try:
exec(open(Ortho4XP_dir+dir_sep+'Help.py').read(),globals())
except:
pass
for item in configvars:
try:
globals()[item]=configvars_defaults[item]
except:
print("I could not set the variable",item,". Perhaps was there a typo ?")
# These are not put in the interface
build_dir = "default"
tricky_provider_hack= 70000 # The minimum size a wms2048 image should be to be accepted (trying to avoid missed cached)
water_overlay = True
wms_timeout = 60
pools_max_points = 65536 # do not change this !
shutdown_timer = 60 # Time in seconds to close program / shutdown computer after completition
shutd_msg_interval = 15 # Shutdown message display interval
raster_resolution = 10000 # Image size for the raster of the sniffed landclass terrain, not yet used
# Will be used as global variables
download_to_do_list=[]
convert_to_do_list=[]
busy_slots_mont=0
busy_slots_conv=0
if 'dar' in sys.platform:
dir_sep = '/'
Triangle4XP_cmd = Ortho4XP_dir+"/Utils/Triangle4XP.app "
copy_cmd = "cp "
delete_cmd = "rm "
rename_cmd = "mv "
unzip_cmd = "7z "
convert_cmd = "convert "
convert_cmd_bis = Ortho4XP_dir+dir_sep+"Utils"+dir_sep+"nvcompress"+dir_sep+"nvcompress-osx -bc1 -fast "
gimp_cmd = "gimp "
devnull_rdir = " >/dev/null 2>&1"
shutdown_cmd = 'sudo shutdown -h now'
os.system('chmod a+x '+Ortho4XP_dir+dir_sep+'Utils/DSFTool.app')
os.system('chmod a+x '+Ortho4XP_dir+dir_sep+'Utils/Triangle4XP.app')
os.system('chmod a+x '+Ortho4XP_dir+dir_sep+'Utils/nvcompress/nvcompress-osx')
elif 'win' in sys.platform:
dir_sep = '\\'
Triangle4XP_cmd = Ortho4XP_dir+dir_sep+"Utils"+dir_sep+"Triangle4XP.exe "
copy_cmd = "copy "
delete_cmd = "del "
rename_cmd = "move "
unzip_cmd = Ortho4XP_dir+dir_sep+"Utils"+dir_sep+"7z.exe "
if os.path.isfile(Ortho4XP_dir+dir_sep+"Utils"+dir_sep+"convert.exe"):
convert_cmd = Ortho4XP_dir+dir_sep+"Utils"+dir_sep+"convert.exe "
else:
convert_cmd = "convert "
convert_cmd_bis = Ortho4XP_dir+dir_sep+"Utils"+dir_sep+"nvcompress"+dir_sep+"nvcompress.exe -bc1 -fast "
gimp_cmd = "c:\\Program Files\\GIMP 2\\bin\\gimp-console-2.8.exe "
showme_cmd = Ortho4XP_dir+"/Utils/showme.exe "
devnull_rdir = " > nul 2>&1"
shutdown_cmd = 'shutdown /s /f /t 0'
else:
dir_sep = '/'
Triangle4XP_cmd = Ortho4XP_dir+"/Utils/Triangle4XP "
delete_cmd = "rm "
copy_cmd = "cp "
rename_cmd = "mv "
unzip_cmd = "7z "
convert_cmd = "convert "
convert_cmd_bis = "nvcompress -fast -bc1a "
gimp_cmd = "gimp "
devnull_rdir = " >/dev/null 2>&1 "
shutdown_cmd = 'sudo shutdown -h now'
os.system('chmod a+x '+Ortho4XP_dir+dir_sep+'Utils/DSFTool')
os.system('chmod a+x '+Ortho4XP_dir+dir_sep+'Utils/Triangle4XP')
#
# END OF FACTORY DEFAULT VALUES
#
##############################################################################
dico_edge_markers = {'outer':'1','inner':'1','coastline':'2',\
'tileboundary':'3','orthogrid':'3',\
'airport':'4','runway':'5','patch':'6'}
dico_tri_markers = {'water':'1','sea':'2','sea_equiv':'3','altitude':'4'}
try:
exec(open(Ortho4XP_dir+dir_sep+'Carnet_d_adresses.py').read())
except:
print("The file Carnet_d_adresses.py does not follow the syntactic rules.")
time.sleep(5)
sys.exit()
try:
exec(open(Ortho4XP_dir+dir_sep+'APL_scripts.py').read())
except:
pass
##############################################################################
# Minimalist error messages. #
##############################################################################
##############################################################################
def usage(reason,do_i_quit=True):
if reason=='config':
print("The file Ortho4XP.cfg was not found or does not follow the "+\
"syntactic rules.")
elif reason=='command_line':
print("The command line does not follow the syntactic rules.")
elif reason=='osm_tags':
print("I had a problem downloadings data from Openstreetmap.\n"+\
"Your connection may be unavailable or the Overpass server\n"+\
"may be unreachable.")
elif reason=='dem_files':
print("!!!I could not fin the elevation data file, or it was broken.\n!!!I go on with all zero elevation (perhaps a tile full of sea ?)")
elif reason=='adresses':
print("The file Carnet_d_adresses.py does not follow the syntactic"+\
" rules.")
elif reason=='crash':
print("The mesh algorithm Triangle4XP has encountered a problem and"+\
" had to stop.")
elif reason=='inprogress':
print("This functionality is not yet supported.")
if do_i_quit==True:
sys.exit()
return
##############################################################################
#######################
#
# For future use
#
def build_landclass_poly_file(lat0,lon0,build_dir,file_to_sniff,dem_alternative=False):
re_encode_dsf(lat0,lon0,build_dir,file_to_sniff,True,dem_alternative,True)
return
# Don't go further now !
t1=time.time()
if not os.path.exists(build_dir):
os.makedirs(build_dir)
strlat='{:+.0f}'.format(lat0).zfill(3)
strlon='{:+.0f}'.format(lon0).zfill(4)
poly_file = build_dir+dir_sep+'Data'+strlat+strlon+'.poly'
airport_file = build_dir+dir_sep+'Data'+strlat+strlon+'.apt'
patch_dir = Ortho4XP_dir+dir_sep+'Patches'+dir_sep+strlat+strlon
dico_nodes={}
dico_edges={}
dico_edges_water={}
water_seeds=[]
sea_seeds=[]
sea_equiv_seeds=[]
flat_airport_seeds=[]
flat_patch_seeds=[]
sloped_patch_seeds=[]
alt_seeds=[]
zone_list=[] # !!!!!!!!!!!!!!!!!!!!!!!!
poly_list=[region[0] for region in zone_list]
bbox_list=[compute_bbox(poly) for poly in poly_list]
print("-> Analyzing patch data ")
include_patch_data(lat0,lon0,patch_dir,dico_nodes,dico_edges,flat_patch_seeds,sloped_patch_seeds,alt_seeds,poly_list,bbox_list)
print(poly_list)
print("-> Analyzing input landclass mesh ")
(tri_list,ter_list,tri_list_kept,hole_seeds,textures)=read_dsf(lat0,lon0,build_dir,file_to_sniff,poly_list,bbox_list,dem_alternative)
return textures
#hole_seeds=[]
print("tri_list :",len(tri_list))
print("ter_list :",len(ter_list))
#for key in tri_list_kept:
# try:
# print(key)
# print(len(tri_list_kept[key]))
# print(' ')
# except:
# pass
#print("-> Computing splitting edges and building raster landclass map ")
#tri_list=split_mesh_according_to_poly_list(lat0,lon0,triangles,poly_list,bbox_list)
print("-> Building poly file for Triangle4XP ")
ttest=time.time()
for tri,ter in zip(tri_list,ter_list):
node1=['{:.6f}'.format(tri[0][1]),'{:.6f}'.format(tri[0][0])] # i.e. [lat,lon]
node2=['{:.6f}'.format(tri[1][1]),'{:.6f}'.format(tri[1][0])] # i.e. [lat,lon]
node3=['{:.6f}'.format(tri[2][1]),'{:.6f}'.format(tri[2][0])] # i.e. [lat,lon]
keep_node(node1,lat0,lon0,dico_nodes,tri[0][2])
keep_node(node2,lat0,lon0,dico_nodes,tri[1][2])
keep_node(node3,lat0,lon0,dico_nodes,tri[2][2])
keep_edge_unique(node1,node2,'inner',dico_edges)
keep_edge_unique(node2,node3,'inner',dico_edges)
keep_edge_unique(node3,node1,'inner',dico_edges)
if ter==0: # water
keep_edge(node1,node2,'inner',dico_edges_water)
keep_edge(node2,node3,'inner',dico_edges_water)
keep_edge(node3,node1,'inner',dico_edges_water)
#if tri[3]==0: # terrain_Water
# tri_bary=[(tri[0][0]+tri[1][0]+tri[2][0])/3,(tri[0][1]+tri[1][1]+tri[2][1])/3]
# water_seeds.append(tri_bary)
for edge in dico_edges_water:
node0=edge.split('|')[0].split('_')
node1=edge.split('|')[1].split('_')
keep_edge(node0,node1,'inner',dico_edges)
print("time to keep edges : ",time.time()-ttest)
#if zone_list != []:
# print("-> Adding of edges related to the orthophoto grid and computation of\n"
# " their intersections with OSM edges,")
# dico_edges=cut_edges_with_grid(lat0,lon0,dico_nodes,dico_edges)
print(" Removal of obsolete edges,")
dico_edges_tmp={}
for edge in dico_edges:
[initpt,endpt]=edge.split('|')
if initpt != endpt:
dico_edges_tmp[edge]=dico_edges[edge]
else:
print("one removed edge : "+str(initpt))
dico_edges=dico_edges_tmp
print(" Removal of obsolete nodes,")
final_nodes={}
for edge in dico_edges:
#print(edge)
[initpt,endpt]=edge.split('|')
final_nodes[initpt]=dico_nodes[initpt]
final_nodes[endpt]=dico_nodes[endpt]
dico_nodes=final_nodes
print("-> Transcription of the data to the file "+poly_file)
write_poly_file(poly_file,airport_file,lat0,lon0,dico_nodes,dico_edges, water_seeds,sea_seeds,sea_equiv_seeds,\
flat_airport_seeds,flat_patch_seeds,sloped_patch_seeds,alt_seeds,hole_seeds)
print('\nCompleted in '+str('{:.2f}'.format(time.time()-t1))+\
'sec.')
print('_____________________________________________________________'+\
'____________________________________')
return
##############################################################################
# Construction du fichier .poly décrivant toutes les données vectorielles
# a intégrer au maillage (frontières sol/eau et aéroports).
##############################################################################
def build_poly_file(lat0,lon0,option,build_dir,airport_website=default_website):
t1=time.time()
if not os.path.exists(build_dir):
os.makedirs(build_dir)
strlat='{:+.0f}'.format(lat0).zfill(3)
strlon='{:+.0f}'.format(lon0).zfill(4)
if not os.path.exists(Ortho4XP_dir+dir_sep+"OSM_data"+dir_sep+strlat+strlon):
os.makedirs(Ortho4XP_dir+dir_sep+"OSM_data"+dir_sep+strlat+strlon)
poly_file = build_dir+dir_sep+'Data'+strlat+strlon+'.poly'
airport_file = build_dir+dir_sep+'Data'+strlat+strlon+'.apt'
patch_dir = Ortho4XP_dir+dir_sep+'Patches'+dir_sep+strlat+strlon
dico_nodes={}
dico_edges={}
# seeds will serve to populate corresponding regions with the appropriate
# marker
water_seeds=[]
sea_seeds=[]
sea_equiv_seeds=[]
flat_airport_seeds=[]
flat_patch_seeds=[]
sloped_patch_seeds=[]
alt_seeds=[]
init_nodes=0
tags=[]
if option==2: # Orthophoto only for inland water
print("-> Downloading airport and water/ground boundary data on Openstreetmap")
tags.append('way["aeroway"="aerodrome"]')
tags.append('rel["aeroway"="aerodrome"]')
tags.append('way["aeroway"="heliport"]')
tags.append('way["natural"="coastline"]')
else: # Mixed
print("-> Downloading airport and water/ground boundary data from Openstreetmap :")
tags.append('way["aeroway"="aerodrome"]')
tags.append('rel["aeroway"="aerodrome"]')
tags.append('way["aeroway"="heliport"]')
tags.append('way["natural"="water"]["tidal"!="yes"]')
tags.append('rel["natural"="water"]["tidal"!="yes"]')
tags.append('way["waterway"="riverbank"]')
tags.append('rel["waterway"="riverbank"]')
tags.append('way["natural"="coastline"]')
tags.append('way["waterway"="dock"]')
try:
application.red_flag.set(0)
except:
pass
for tag in tags:
try:
if application.red_flag.get()==1:
print("\nOSM download process interrupted.")
print('_____________________________________________________________'+\
'____________________________________')
return
except:
pass
subtags=tag.split('"')
osm_filename=Ortho4XP_dir+dir_sep+"OSM_data"+dir_sep+strlat+strlon+dir_sep+strlat+strlon+'_'+subtags[0][0:-1]+'_'+\
subtags[1]+'_'+subtags[3]+'.osm'
osm_errors_filename=Ortho4XP_dir+dir_sep+"OSM_data"+dir_sep+strlat+strlon+dir_sep+strlat+strlon+'_'+subtags[0][0:-1]+'_'+\
subtags[1]+'_'+subtags[3]+'_detected_errors.txt'
if not os.path.isfile(osm_filename):
print(" Obtaining OSM data for "+tag)
s=requests.Session()
osm_download_ok = False
while osm_download_ok != True:
url=overpass_server_list[overpass_server_choice]+"?data=("+tag+"("+str(lat0)+","+str(lon0)+","+str(lat0+1)+","+str(lon0+1)+"););(._;>>;);out meta;"
r=s.get(url)
if r.headers['content-type']=='application/osm3s+xml':
osm_download_ok=True
else:
print(" OSM server was busy, new tentative...")
osmfile=open(osm_filename,'wb')
osmfile.write(r.content)
osmfile.close()
print(" Done.")
else:
print(" Recycling OSM data for "+tag)
if 'way[' in tag:
[dicosmw,dicosmw_name,dicosmw_icao,dicosmw_ele]=osmway_to_dicos(osm_filename)
elif 'rel[' in tag:
[dicosmr,dicosmrinner,dicosmrouter,dicosmr_name,dicosmr_icao,dicosmr_ele]=osmrel_to_dicos(osm_filename,osm_errors_filename)
# we shall treat osm data differently depending on tag
if tag=='way["aeroway"="aerodrome"]' or tag=='way["aeroway"="heliport"]':
sloped_airports_list=[]
if os.path.exists(patch_dir):
for pfilename in os.listdir(patch_dir):
if (pfilename[-10:] == '.patch.osm') or os.path.isdir(patch_dir+dir_sep+pfilename):
sloped_airports_list.append(pfilename[:4])
#print(sloped_airports_list)
for wayid in dicosmw:
way=dicosmw[wayid]
# we only treat closed ways, non closed should not exist in
# osm ways (but a few sometimes do!)
if strcode(way[0]) == strcode(way[-1]):
signed_area=area(way)
if signed_area<0:
side='left'
else:
side='right'
keep_that_one=True
if (wayid in dicosmw_icao) and (wayid in dicosmw_name):
print(" * "+dicosmw_icao[wayid]+" "+dicosmw_name[wayid])
elif (wayid in dicosmw_icao):
print(" * "+dicosmw_icao[wayid])
elif (wayid in dicosmw_name):
print(" * "+dicosmw_name[wayid])
else:
print(" * lat="+str(way[0][0])+" lon="+str(way[0][1]))
if (wayid in dicosmw_ele):
altitude=dicosmw_ele[wayid]
else:
altitude='unknown'
if (wayid in dicosmw_icao):
if (dicosmw_icao[wayid] in sloped_airports_list) or (dicosmw_icao[wayid] in do_not_flatten_these_list):
print(" I will not flatten "+dicosmw_icao[wayid]+ " airport.")
keep_that_one=False
if keep_that_one==True:
keep_way(way,lat0,lon0,1,'airport',dico_nodes,\
dico_edges)
flat_airport_seeds.append([way,\
pick_point_check(way,side,lat0,lon0),altitude])
else:
print("One of the airports within the tile is not correctly closed \n"+\
"on Openstreetmap ! Close to coordinates " + str(way[0]))
elif tag=='rel["aeroway"="aerodrome"]':
sloped_airports_list=[]
if os.path.exists(patch_dir):
for pfilename in os.listdir(patch_dir):
if pfilename[-10:] == '.patch.osm':
sloped_airports_list.append(pfilename[:4])
for relid in dicosmr:
keep_that_one=True
if (relid in dicosmr_icao):
if dicosmr_icao[relid] in sloped_airports_list or (dicosmr_icao[relid] in do_not_flatten_these_list):
keep_that_one=False
if (relid in dicosmr_icao) and (relid in dicosmr_name):
print(" * "+dicosmr_icao[relid]+" "+dicosmr_name[relid])
elif (relid in dicosmr_icao):
print(" * "+dicosmr_icao[relid])
elif relid in dicosmr_name:
print(" * "+dicosmr_name[relid])
if (relid in dicosmr_ele):
altitude=dicosmr_ele[relid]
else:
altitude='unknown'
if keep_that_one==False:
continue
for waypts in dicosmrinner[relid]:
keep_way(waypts,lat0,lon0,1,'airport',dico_nodes,\
dico_edges)
for waypts in dicosmrouter[relid]:
signed_area=area(waypts)
if signed_area<0:
side='left'
else:
side='right'
keep_way(waypts,lat0,lon0,1,'airport',dico_nodes,\
dico_edges)
flat_airport_seeds.append([waypts,\
pick_point_check(waypts,side,lat0,lon0),altitude])
elif 'way["natural"="coastline"]' in tag:
total_sea_seeds=0
for wayid in dicosmw:
way=dicosmw[wayid]
# Openstreetmap ask that sea is to the right of oriented
# coastline. We trust OSM contributors...
if strcode(way[0])!=strcode(way[-1]):
if (lat0>=40 and lat0<=49 and lon0>=-93 and lon0<=-76):
sea_equiv_seeds+=pick_points_safe(way,'right',lat0,lon0)
else:
sea_seeds+=pick_points_safe(way,'right',lat0,lon0)
total_sea_seeds+=1
keep_way(way,lat0,lon0,1,'coastline',dico_nodes,dico_edges)
if total_sea_seeds<=3:
for wayid in dicosmw:
way=dicosmw[wayid]
if strcode(way[0])==strcode(way[-1]):
if (lat0>=40 and lat0<=49 and lon0>=-93 and lon0<=-76):
sea_equiv_seeds+=pick_points_safe(way,'right',lat0,lon0)
else:
sea_seeds+=pick_points_safe(way,'right',lat0,lon0)
elif ('way["natural"="water"]' in tag) or ('way["waterway"="riverbank"]' in tag) or ('way["waterway"="dock"]' in tag) :
efile=open(osm_errors_filename,'w')
osm_errors_found=False
for wayid in dicosmw:
#print(wayid)
way=dicosmw[wayid]
if strcode(way[0]) != strcode(way[-1]):
osm_errors_found=True
efile.write("Way id="+str(wayid)+" was not treated because it is not closed.\n")
continue
if touches_region(way,lat0,lat0+1,lon0,lon0+1):
signed_area=area(way)
# Keep only sufficiently large water pieces. 1 deg^2 is
# very roughly equal to 10000 km^2
if abs(signed_area) >= min_area/10000.0:
if signed_area<0:
side='left'
else:
side='right'
keep_way(way,lat0,lon0,1,'outer',dico_nodes,\
dico_edges)
points_checked=pick_points_safe(way,side,lat0,lon0,check=True)
#print(points_checked)
sea_way=False
try:
if dicosmw_name[wayid] in sea_equiv:
print(" Sea_equiv found :",dicosmw_name[wayid])
sea_way=True
except:
pass
if sea_way!=True:
water_seeds+=points_checked
else:
sea_equiv_seeds+=points_checked
efile.close()
if osm_errors_found:
print(" !!!Some OSM errors were detected!!!\n They are listed in "+str(osm_errors_filename))
else:
os.remove(osm_errors_filename)
elif 'rel[' in tag:
for relid in dicosmr:
sea_rel=False
try:
if dicosmr_name[relid] in sea_equiv:
print(" Sea_equiv found :",dicosmr_name[relid])
sea_rel=True
except:
pass
for waypts in dicosmrinner[relid]:
keep_way(waypts,lat0,lon0,1,'inner',dico_nodes,\
dico_edges)
for waypts in dicosmrouter[relid]:
signed_area=area(waypts)
if abs(signed_area) >= min_area/10000.0:
if signed_area<0:
side='left'
else:
side='right'
keep_way(waypts,lat0,lon0,1,'outer',dico_nodes,\
dico_edges)
points_checked=pick_points_safe(waypts,side,lat0,lon0,check=True)
if sea_rel!=True:
water_seeds+=points_checked
else:
sea_equiv_seeds+=points_checked
treated_nodes=len(dico_nodes)-init_nodes
init_nodes=len(dico_nodes)
if treated_nodes>0:
strtmp=str(treated_nodes)+" new nodes."
else:
strtmp="no new node."
print(" -> process of the associated data completed : "+strtmp)
print("-> Cutting off of too long edges,")
dico_edges_tmp={}
for edge in dico_edges:
[initpt,endpt]=edge.split('|')
[xi,yi]=xycoords(initpt,lat0,lon0)
[xf,yf]=xycoords(endpt,lat0,lon0)
#length=sqrt((xi-xf)*(xi-xf)+(yi-yf)*(yi-yf))
length=abs(xi-xf)+abs(yi-yf)
pieces=ceil(length*1000)
if pieces == 1:
dico_edges_tmp[edge]=dico_edges[edge]
else:
coordlist=[]
for k in range(1,pieces):
xk=((pieces-k)/pieces)*xi+(k/pieces)*xf
yk=((pieces-k)/pieces)*yi+(k/pieces)*yf
coordlist.append([xk,yk])
keep_node_xy(xk,yk,lat0,lon0,dico_nodes)
keep_edge_str_tmp(initpt,strxy(coordlist[0][0],coordlist[0][1],\
lat0,lon0),dico_edges[edge],dico_edges_tmp)
for k in range(1,pieces-1):
keep_edge_str_tmp(strxy(coordlist[k-1][0],coordlist[k-1][1],\
lat0,lon0),strxy(coordlist[k][0],coordlist[k][1],\
lat0,lon0),dico_edges[edge],dico_edges_tmp)
keep_edge_str_tmp(strxy(coordlist[pieces-2][0],\
coordlist[pieces-2][1],lat0,lon0),endpt,dico_edges[edge],\
dico_edges_tmp)
dico_edges=dico_edges_tmp
print("-> Adding patch data for the mesh, ")
include_patch_data(lat0,lon0,patch_dir,dico_nodes,dico_edges,flat_patch_seeds,sloped_patch_seeds,alt_seeds)
if zone_list and insert_custom_zoom_in_mesh:
print("-> Adding of edges related to the custom zoomlevel and computation of\n"
" their intersections with OSM edges,")
dico_edges=cut_edges_with_zone(lat0,lon0,dico_nodes,dico_edges,zone_list)
print("-> Adding of edges related to the orthophoto grid and computation of\n"
" their intersections with OSM edges,")
dico_edges=cut_edges_with_grid(lat0,lon0,dico_nodes,dico_edges)
print(" Removal of obsolete edges,")
dico_edges_tmp={}
for edge in dico_edges:
[initpt,endpt]=edge.split('|')
if initpt != endpt:
dico_edges_tmp[edge]=dico_edges[edge]
#else:
# print("one removed edge : "+str(initpt))
dico_edges=dico_edges_tmp
print(" Removal of obsolete nodes,")
final_nodes={}
for edge in dico_edges:
#print(edge)
[initpt,endpt]=edge.split('|')
final_nodes[initpt]=dico_nodes[initpt]
final_nodes[endpt]=dico_nodes[endpt]
dico_nodes=final_nodes
print("-> Transcription of the updated data to the file "+poly_file)
write_poly_file(poly_file,airport_file,lat0,lon0,dico_nodes,dico_edges, water_seeds,sea_seeds,sea_equiv_seeds,\
flat_airport_seeds,flat_patch_seeds,sloped_patch_seeds,alt_seeds)
# adding airports high zl zone to zone_list if needed
try:
if cover_airports_with_highres==True: #and (zone_list==[] or (zone_list[-1][0][0]==zone_list[-1][0][2]==zone_list[-1][0][8])): # i.e. airports claimed but not already in zone_list (since write_cfg puts them there)
filed_some=False
with open(build_dir+dir_sep+"Data"+strlat+strlon+'.apt') as fapt:
for line in fapt:
if 'Airport' in line: # or 'patch' in line:
nbr_points=int(line.split()[3])
skip=fapt.readline()
latmin=90
latmax=-90
lonmin=180
lonmax=-180
for k in range(0,nbr_points):
[latp,lonp]=fapt.readline().split()
latp=float(latp)
lonp=float(lonp)
latmin = latp if latp<latmin else latmin
latmax = latp if latp>latmax else latmax
lonmin = lonp if lonp<lonmin else lonmin
lonmax = lonp if lonp>lonmax else lonmax
half_meridian=pi*6378137
texsize=int(2*half_meridian/2**(cover_zl-4)*cos(latmin*pi/180))
e=int(round(cover_extent*1000.0/texsize))
[a,b]=wgs84_to_texture(latmax,lonmin,cover_zl,'BI')
[c,d]=wgs84_to_texture(latmin,lonmax,cover_zl,'BI')
[latmax,lonmin]=gtile_to_wgs84(a-16*e,b-16*e,cover_zl)
[latmin,lonmax]=gtile_to_wgs84(c+16*(1+e),d+16*(1+e),cover_zl)
new_entry=[[latmin,lonmin,latmin,lonmax,latmax,lonmax,latmax,lonmin,latmin,lonmin],str(cover_zl),str(airport_website)]
if new_entry not in zone_list:
zone_list.append(new_entry)
filed_some=True
if filed_some: print('-> Airports and patched areas have been added to the custom zoomlevel list with ZL' +str(cover_zl)+' and approx radius of '+str(cover_extent)+' km.')
except:
pass
print('\nCompleted in '+str('{:.2f}'.format(time.time()-t1))+\
'sec.')
print('_____________________________________________________________'+\
'____________________________________')
return
#############################################################################
#############################################################################
def write_poly_file(poly_file,airport_file,lat0,lon0,dico_nodes,dico_edges,\
water_seeds,sea_seeds,sea_equiv_seeds,flat_airport_seeds,\
flat_patch_seeds, sloped_patch_seeds,alt_seeds,hole_seeds=None):
total_nodes=len(dico_nodes)
f=open(poly_file,'w')
f.write(str(total_nodes)+' 2 1 0\n')
dico_node_pos={}
idx=1
for key in dico_nodes:
dico_node_pos[key]=idx
[x,y]=xycoords(key,lat0,lon0)
f.write(str(idx)+' '+str(x)+' '+\
str(y)+' '+str(dico_nodes[key])+'\n')
idx+=1
f.write('\n')
idx=1
total_edges=len(dico_edges)
f.write(str(total_edges)+' 1\n')
for edge in dico_edges:
[code1,code2]=edge.split('|')
idx1=dico_node_pos[code1]
idx2=dico_node_pos[code2]
f.write(str(idx)+' '+str(idx1)+' '+str(idx2)+' '+\
dico_edge_markers[dico_edges[edge]]+'\n')
idx+=1
if hole_seeds==None or len(hole_seeds)==0:
f.write('\n0\n')
else:
nbr_hole_seeds=len(hole_seeds)
f.write('\n'+str(nbr_hole_seeds)+'\n')
idx=1
for seed in hole_seeds:
f.write(str(idx)+' '+str(seed[1]-lon0)+' '+str(seed[0]-lat0)+'\n')
idx+=1
total_seeds=len(water_seeds)+len(sea_seeds)+len(sea_equiv_seeds)+\
len(flat_airport_seeds)+len(flat_patch_seeds)+\
len(sloped_patch_seeds)+len(alt_seeds)
if total_seeds==0:
print("-> No OSM data for this tile, loading altitude data to decide between sea or land.")
[alt_dem,ndem]=load_altitude_matrix(lat0,lon0)
land =(alt_dem.max()-alt_dem.min())>=20
if land: #sea_texture_params==[]:
print(" Decided for a land tile.")
water_seeds.append([1000,1000])
else:
print(" Decided for a sea tile.")
sea_seeds.append([lon0+0.5,lat0+0.5])
total_seeds=1
f.write('\n'+str(total_seeds)+' 1\n')
idx=1
for seed in water_seeds:
f.write(str(idx)+' '+str(seed[0]-lon0)+' '+str(seed[1]-lat0)+' '+\
dico_tri_markers['water']+'\n')
idx+=1
for seed in sea_seeds:
f.write(str(idx)+' '+str(seed[0]-lon0)+' '+str(seed[1]-lat0)+' '+\
dico_tri_markers['sea']+'\n')
idx+=1
for seed in sea_equiv_seeds:
f.write(str(idx)+' '+str(seed[0]-lon0)+' '+str(seed[1]-lat0)+' '+\
dico_tri_markers['sea_equiv']+'\n')
idx+=1
for seed in alt_seeds:
f.write(str(idx)+' '+str(seed[0]-lon0)+' '+str(seed[1]-lat0)+' '+\
dico_tri_markers['altitude']+'\n')
idx+=1
apt_idx=100
for seed in flat_airport_seeds:
f.write(str(idx)+' '+str(seed[1][0]-lon0)+' '+\
str(seed[1][1]-lat0)+' '+str(apt_idx)+'\n')
apt_idx+=1
idx+=1
fp_idx=1000
for seed in flat_patch_seeds:
f.write(str(idx)+' '+str(seed[0][0])+' '+\
str(seed[0][1])+' '+str(fp_idx)+'\n')
fp_idx+=1
idx+=1
sp_idx=10000
for seed in sloped_patch_seeds:
f.write(str(idx)+' '+str(seed[0][0])+' '+str(seed[0][1])+' '+\
str(sp_idx)+'\n')
sp_idx+=1
idx+=1
print(" Remain " + str(len(dico_edges))+\
" edges in total.")
f.close()
f=open(airport_file,"w")
apt_idx = 100
fp_idx = 1000
sp_idx = 10000
for seed in flat_airport_seeds:
f.write("Airport "+str(apt_idx)+" : "+str(len(seed[0]))+\
" nodes.\n")
f.write("Elevation "+str(seed[2])+'\n')
for node in seed[0]:
f.write(str(float(node[0]))+" "+str(float(node[1]))+"\n")
f.write("\n")
apt_idx+=1
f.write('\n')
for seed in flat_patch_seeds:
f.write("Flat_patch "+str(fp_idx)+" : "+str(len(seed[2]))+"\n")
f.write("Elevation "+str(seed[1])+'\n')
for node in seed[2]:
[slat,slon]=node.split('_')
f.write(slat+" "+slon+"\n")
f.write("\n")
fp_idx+=1
f.write('\n')
for seed in sloped_patch_seeds:
f.write("Sloped_patch "+str(sp_idx)+" : "+str(seed[1])+" "+\
str(seed[2])+" "+str(seed[3])+" "+str(seed[4])+" "+\
str(seed[5])+" "+str(seed[6])+" "+str(seed[7])+" "+\
str(seed[8])+" "+str(seed[9])+"\n")
sp_idx+=1
f.close()
return
##############################################################################
##############################################################################
def osmway_to_dicos(osm_filename):
pfile=open(osm_filename,'r',encoding="utf-8")
dicosmn={}
dicosmw={}
dicosmw_name={}
dicosmw_icao={}
dicosmw_ele={}
finished_with_file=False
in_way=False
first_line=pfile.readline()
if "'" in first_line:
separator="'"
else:
separator='"'
while not finished_with_file==True:
items=pfile.readline().split(separator)
if '<node id=' in items[0]:
id=items[1]
for j in range(0,len(items)):
if items[j]==' lat=':
slat=items[j+1]
elif items[j]==' lon=':
slon=items[j+1]
dicosmn[id]=[slat,slon]
elif '<way id=' in items[0]:
in_way=True
wayid=items[1]
dicosmw[wayid]=[]
elif '<nd ref=' in items[0]:
dicosmw[wayid].append(dicosmn[items[1]])
elif '<tag k=' in items[0] and in_way and items[1]=='name':
dicosmw_name[wayid]=items[3]
elif '<tag k=' in items[0] and in_way and items[1]=='icao':
dicosmw_icao[wayid]=items[3]
elif '<tag k=' in items[0] and in_way and items[1]=='ele':
dicosmw_ele[wayid]=items[3]
elif '</osm>' in items[0]:
finished_with_file=True
pfile.close()
# (Thanks Tony Wroblewski) Remove ways with no nodes (Can happen when editing in JOSM or from bad data)
for wayid in list(dicosmw.keys()):
if not dicosmw[wayid]:
del dicosmw[wayid]
print(" A total of "+str(len(dicosmn))+" node(s) and "+str(len(dicosmw))+" way(s).")
return [dicosmw,dicosmw_name,dicosmw_icao,dicosmw_ele]
##############################################################################
##############################################################################
def osmrel_to_dicos(osm_filename,osm_errors_filename):
pfile=open(osm_filename,'r',encoding="utf-8")
efile=open(osm_errors_filename,'w')
osm_errors_found=False
dicosmn={}
dicosmw={}
dicosmr={}
dicosmrinner={}
dicosmrouter={}
dicosmr_name={}
dicosmr_icao={}
dicosmr_ele={}
dicoendpt={}
finished_with_file=False
in_rel=False
first_line=pfile.readline()
if "'" in first_line:
separator="'"
else:
separator='"'
while not finished_with_file==True:
items=pfile.readline().split(separator)
if '<node id=' in items[0]:
id=items[1]
for j in range(0,len(items)):
if items[j]==' lat=':
slat=items[j+1]
elif items[j]==' lon=':
slon=items[j+1]
dicosmn[id]=[slat,slon]
elif '<way id=' in items[0]:
wayid=items[1]
dicosmw[wayid]=[]
elif '<nd ref=' in items[0]:
dicosmw[wayid].append(items[1])
elif '<relation id=' in items[0]:
relid=items[1]
in_rel=True
dicosmr[relid]=[]
dicosmrinner[relid]=[]
dicosmrouter[relid]=[]
dicoendpt={}
elif '<member type=' in items[0]:
if items[1]!='way':
efile.write("Relation id="+str(relid)+" contains member "+items[1]+" which was not treated because it is not a way.\n")
osm_errors_found=True
continue
role=items[5]
if role=='inner':
waytmp=[]
for nodeid in dicosmw[items[3]]:
waytmp.append(dicosmn[nodeid])
dicosmrinner[relid].append(waytmp)
elif role=='outer':
endpt1=dicosmw[items[3]][0]
endpt2=dicosmw[items[3]][-1]
if endpt1==endpt2:
waytmp=[]
for nodeid in dicosmw[items[3]]:
waytmp.append(dicosmn[nodeid])
dicosmrouter[relid].append(waytmp)
else:
if endpt1 in dicoendpt: