Mer

Rasterfliser fra vektordata med GDAL - hvordan unngå å aliassere gjenstander?


Prequesites

Jeg har vektordata i en ShapefileN47E007_CONTOURS.shpprodusert med

wget -O N47E007.zip http://e4ftl01.cr.usgs.gov/SRTM/SRTMGL1.003/2000.02.11/N47E007.SRTMGL1.hgt.zip pakke ut N47E007.zip # Konverter fra WGS 84 til Web Mercator og fra HGT til GeoTIFF gdalwarp -s_srs EPSG: 4326 N47E007.hgt  -t_srs EPSG: 3785 N47E007_WebMercator.tif  -co "COMPRESS = LZW" -av GTiff -r bilineær # Beregn konturvektordata gdal_contour -a elev -i 20 N47E007_ONT .shp

(Konturlinjer basert på en lapp (N47E007) av NASA SRTM høydedata.)

Prosess

Fra dette genererer jeg TMS / WMTS-fliser som skal brukes av brosjyren, ved først å rastere vektordataene fra Shapefile medgdal_rasterizeog deretter flislegge den medgdal2tiles.py. Jeg gjentar begge trinnene for hvert zoomnivå, slik at jeg kan rasterisere til en større raster for høyere zoomnivåer, slik at linjene får samme tykkelse på skjermen på hver enkelt:

#! / bin / sh for z i $ (seq 1 15) do gdal_rasterize  --config GDAL_CACHEMAX 1024  -ts $ (calc "2851 << ($ {z} - 12)")  $ (calc "4220 < <($ {z} - 12) ")  -ot byte -of GTiff -co COMPRESS = DEFLATE -co PREDICTOR = 2 -co ZLEVEL = 9  -ot Byte -burn 0 -a_nodata 255  -l N47E007_CONTOURS N47E007_CONTOURS.shp N47E007_CONTOURS_rasterized_z $ {z} .tif gdal2tiles.py  -r nær -z $ {z} -a 255255255  -e N47E007_CONTOURS_rasterized_z $ {z} .tif  N47E007_contours_tiled gjort

Tilpasser oppløsningen til zoomnivået

(Se også svaret mitt på Problem med konturlinjetykkelse i større zoomnivåer.)

beregne "2851 << ($ {z} - 12)"ogberegne "4220 << ($ {z} - 12)"beregne størrelsen på bildet du vil raste til (N47E007_CONTOURS_rasterized_z $ {z} .tifavhengig av zoomnivået$ z.

<<er den binære venstre-skiftoperatøren, så dette oversettes til 2851⋅2 z-12 og 4220⋅2 z-12. Dette betyr at

  • for zoomnivå 12, det mellomliggende rasterbildetN47E007_CONTOURS_rasterized_z12.tifvil være 2851 × 4220 piksler,
  • for zoomnivå 13 (N47E007_CONTOURS_rasterized_z13.tif) vil det være dobbelt så mye i hver retning (5702 × 8440) og
  • for zoomnivå 11 (N47E007_CONTOURS_rasterized_z11.tif) halvparten av zoomnivå 12, dvs. 1425 × 2110.

Presis parametrisering

Tallene2851og4220(størrelsen for zoomnivå 12) ble gjettet eller målt på skjermen eller på annen måte avledet på en upresis måte.

Problem

De resulterende flisene viser noen gjenstander som ikke er synlige i mellomliggende (til) rasteriserte bilder:

  • doble bredde segmenter av linjer,
    f.eks.13/4270 / 5331.png "> vs.

  • hull i linjer,
    f.eks. ved flisekanter mellom

    14/8536 / 10665.png ">

    tilsvarende utdrag fraN47E007_CONTOURS_rasterized_z14.tiftil sammenligning:

Hvis jeg hopper over-r nær(så rasterisering bruker standard resampling-modus 'gjennomsnitt') Jeg får i stedet noen utvaskede piksler som får linjene til å se uskarpe ut.

Teorien min

Koordinatsystemene samsvarer allerede (på grunn avgdalwarptrinn før du beregner konturene), så det er ikke nødvendig å kreve generell transformasjon / omprojeksjon. Men fordi størrelsen som skal rasteriseres til bare tilnærmer størrelsen på skjermen til det geografiske området som dekkes avN47E007dataoppdatering, må noen omskalering fortsatt skje igdal2tiles.py, for å produsere sømløse fliser med riktig plassering og størrelse.

Med-r nærdette forårsaker en slags makroskopisk Moiré-effekt, fordi pikselnettene til flisene og den mellomliggende rasteren nesten bare stemmer overens. (I tillegggdal2tiles.pyser kanskje ikke ut over fliser. Når vi ser etter nærmeste piksel i kildebildet, forårsaker hullene.) Withgjennomsnitt(standard for-r), påvirker de samme årsakene det utvaskede utseendet fordi verdier fra mer enn en kildepiksel blir 'blandet' for å få verdien for en gitt flisepiksel.

Dermed, hvis jeg visste den riktige (zoomnivåavhengige) skjermstørrelsen til det geografiske området til datapatchen og vil bruke det til rasterisering, bør artefaktene forsvinne.

Spørsmål

  • Er teorien min riktig?
  • I så fall hvordan kan jeg bestemme nøyaktig-tsparametere forgdal_rasterize(eller alternativt-trparametere) slik at ingen skalering (eller en omskalering uten effekt) vil bli gjort avgdal2tiles.py?

Er teorien min riktig?

Det virker som jeg får feilfrie resultater med metoden beskrevet nedenfor.

I så fall hvordan kan jeg bestemme nøyaktig-tsparametere forgdal_rasterize(eller alternativt-trparametere) slik at ingen skalering (eller en omskalering uten effekt) vil bli gjort avgdal2tiles.py?

Samtidig somgdal2tiles.pyer ment å brukes som en kommando, kan vi laste den som en pythonmodul og bruke de samme funksjonene den bruker selv for størrelsesberegningene. Dermed har jeg laget meg dette lille python-skriptet:

fra __future__ import utskriftsfunksjon import sys fra gdal2tiles import GlobalMercator def main (): lat, lon = [float (arg) for arg i sys.argv [1: 3]] zoom = int (sys.argv [3]) px, py = lat_lon_to_pixels (lat, lon, zoom) print (px, py) def lat_lon_to_pixels (lat, lon, zoom): gm = GlobalMercator () mx, my = gm.LatLonToMeters (lat, lon) retur tuple (int (rund (p )) for p i gm.MetersToPixels (mx, min, zoom)) hvis __name__ == '__main__': main ()

Jeg bruker deretter WGS 84-koordinatene fragdalinfo N47E007.hgt(trekker dem ut med regulære uttrykk/ ^ s * Øvre høyre s * ( s * (? d + .? d *), s * (? d + .? d *) )og/ ^ s * Nederst til venstre s * ( s * (? d + .? d *), s * (? d + .? d *) ) /) å ringelat_lon_to_pixelsfra python-skriptet ovenfor for det nordøstlige hjørnet og det sørvestlige hjørnet og beregne forskjellene, som jeg deretter bruker som argument for-tsmulighet forgdal_rasterize. (Dette fungerer bare fordi denne HGT-filen har WGS 84 som referansesystem for koordinater, og WGS 84 og WebMercator er aksejustert til hverandre.)


Se videoen: Hoe verwerk je ECCOgravel bij het aanleggen van een grindpad of een oprit van grind? TuinVisie (Oktober 2021).