Jak ponownie rzutować raster z 0 360 na -180 180 za pomocą cięcia południka 180

31

Mam obraz rastrowy geotiff, który ma układ współrzędnych o długości od 0 do 360. Poziomy środek obrazu to 180 długości geograficznej. Zobacz zdjęcie poniżej:

wprowadź opis zdjęcia tutaj

Chcę go przekształcić w EPSG: 4326 SRS o zakresie długości -180 180. Chcę, aby środek obrazu znajdował się na południku Greenwich (0). Myślę, że ten pan jest bardzo szeroko stosowany. Oczekuję, że wynik wygląda następująco:

wprowadź opis zdjęcia tutaj

Więc używam polecenia gdalwarp, aby ponownie wykonać rzut:

gdalwarp -s_srs '+proj=latlong +datum=WGS84 +pm=180dW' -t_srs EPSG:4326 test_col.tif test_4326.tif

Ale dostaję tylko tiff z większymi wymiarami (więcej pikseli) i metadanymi EPSG: 4326. Sam obraz wygląda tak samo, jak początkowy. Ale spodziewam się, że zamieni półkule.

Pytanie brzmi - w jaki sposób mogę uzyskać obraz ściśle -180 180 EPSG: 4326 ze środkiem na długości 0?

To jest gdalinfo mojego początkowego pliku:

Origin = (-0.102272598067084,89.946211604095552)
Pixel Size = (0.204545196134167,-0.204423208191126)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=BAND
Corner Coordinates:
Upper Left  (  -0.1022726,  89.9462116) (  0d 6' 8.18"W, 89d56'46.36"N)
Lower Left  (  -0.1022726, -89.9462116) (  0d 6' 8.18"W, 89d56'46.36"S)
Upper Right (     359.897,      89.946) (359d53'50.18"E, 89d56'46.36"N)
Lower Right (     359.897,     -89.946) (359d53'50.18"E, 89d56'46.36"S)
Center      ( 179.8975000,  -0.0000000) (179d53'51.00"E,  0d 0' 0.00"S)

To jest gdalinfo po gdalwarp

Origin = (-180.102727401932952,89.946211604095552)
Pixel Size = (0.091397622896436,-0.091420837939082)
Metadata:
  AREA_OR_POINT=Area
Image Structure Metadata:
  INTERLEAVE=PIXEL
Corner Coordinates:
Upper Left  (-180.1027274,  89.9462116) (180d 6' 9.82"W, 89d56'46.36"N)
Lower Left  (-180.1027274, -89.9699975) (180d 6' 9.82"W, 89d58'11.99"S)
Upper Right ( 179.8211116,  89.9462116) (179d49'16.00"E, 89d56'46.36"N)
Lower Right ( 179.8211116, -89.9699975) (179d49'16.00"E, 89d58'11.99"S)
Center      (  -0.1408079,  -0.0118929) (  0d 8'26.91"W,  0d 0'42.81"S)
nextstopsun
źródło
O innej rozdzielczości, czy próbowałeś dodać -tr xres yresflagę?
nickves

Odpowiedzi:

21

Można jawnie ustawić wyjściowy zakres współrzędnych za pomocą opcji zasięgu docelowego na gdalwarp (tj. „-Te -180 -90 180 90”), ale można również użyć opcji konfiguracji CENTER_LONG, aby wymusić ponowne zawijanie wokół nowej środkowej długości geograficznej. Coś takiego:

  gdalwarp -t_srs WGS84 ~/0_360.tif 180.tif  -wo SOURCE_EXTRA=1000 \
           --config CENTER_LONG 0

Zwróć także uwagę na opcję wypaczenia „SOURCE_EXTRA = 1000”. Podczas przepakowywania obliczanie prostokąta źródłowego będzie mylone wokół przerwania długości geograficznej i utraci trochę zdjęć. Ta opcja mówi, że pociągnij za trochę. Bez tego zobaczysz lukę danych w pobliżu głównego południka.

PS. ustawianie głównego południka na 180dW tak jak ty nie jest dobrym pomysłem IMHO.

Frank Warmerdam
źródło
1
hmm, --config CENTER_LONG 0nic nie robi, wynikiem jest ten sam raster. Coś tu tęsknię? Działa na GDAL w wersji 2.2.3.
jurajb
6

Zasadniczo musisz pociąć raster na dwie części i złożyć je z powrotem z nowym przesunięciem / skalą.

Jest tutaj przykład, jak to zrobić od [-180,180] do [0,360] za pomocą gdal_translate i sterownika VRT: http://trac.osgeo.org/gdal/wiki/UserDocs/RasterProcTutorial

Zeskanuj w dół do „5-minutowego samouczka”, a szczegóły znajdują się w „Plikach wirtualnych”. Powinno być wystarczająco proste, aby zmodyfikować przykład, aby pasował.

mdsumner
źródło
2

Można to zrobić w R za pomocą jednego wiersza kodu, używając rotatefunkcji z rasterpakietem.

library(raster)
your_raster <- raster("path/to/raster.tif")
rotated_raster <- rotate(your_raster)
SoilSciGuy
źródło
1

Jeśli chcesz tylko wyświetlić raster w QGIS, możesz ustawić niestandardową projekcję za pomocą parametru + lon_wrap = 180.

Rozumiem, że domyślnie proj4 otacza szerokości geograficzne od 0 -> 360 do -180 -> 180. + lon_wrap = 180 skutecznie anuluje to zawijanie i wyświetla szerokości geograficzne między 180 a 360 na półkuli zachodniej.

Opcja + powyżej powinna całkowicie wyłączyć zawijanie, ale - przynajmniej w moim przypadku - raster nie wyświetlał się poprawnie, gdy użyto tej opcji.

Aby uzyskać więcej informacji, zobacz http://proj4.org/parameters.html#lon-wrap-over-longitude-wrapping .


źródło
0

Oto funkcja, którą zbudowałem w celu przerzucenia pojedynczej przyciemnionej tablicy wartości siatki za pomocą javascript od 0-360 do -180-180. Mam nadzieję, że może komuś pomóc.

  let xstart = 180 / xres //xres is the number of values per 1 degree
  for (let y = 0; y < data.height; y++) {
    let index = (y * data.width) + 1,
    start = index + xstart,
    end = index + data.width
    array.splice(index, 0, ...array.splice(start, (end - start)))
  }
maeneak
źródło