5

からncファイルをダウンロードしました 

f=open.ncdf("file.nc")
[1] "file Lfile.nc has  2 dimensions:"
[1] "Longitude   Size: 1440"
[1] "Latitude   Size: 720"
[1] "------------------------"
[1] "file filr.nc has   8 variables:"
[1] "short ts[Latitude,Longitude]  Longname:Skin Temperature (2mm) Missval:NA"

次に、変数soil_moisture_cを操作したいと思いました

A = get.var.ncdf(nc=f,varid="soil_moisture_c",verbose=TRUE)

次に、Aをでプロットしimage(A)ます。下の地図を手に入れました。転置もしましimage(t(a))たが、別の方向に変更されており、本来あるべき姿ではありませんでした。とにかく、何が悪いのかを知るために、私はnetcdfビューアPanoplyを使用し、以下に示すようにマップが正しくプロットされました。

4

3 に答える 3

8

その理由は、使用しているNetCDFインターフェイスが非常に低レベルであり、ディメンション情報なしで変数を読み取るだけであるためです。グリッドの方向は実際には任意であり、座標情報は特定のコンテキストで理解する必要があります。

library(raster) ## requires ncdf package for this file  
d <- raster("LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T185959Z_20040114.nc", varname = "soil_moisture_c")

（私はあなたとは異なるファイルを使用しましたが、同じように機能するはずです）。

ラスターでさえ作業なしではこれを正しく行うことはできませんが、修正は簡単です。

d <-  flip(t(d), direction = "x")

これにより、データが転置され、「x」（経度）が反転し、ジオリファレンスが元のコンテキストから維持されます。

maptoolsからのマップを使用してプロットし、以下を確認します。

plot(d)

library(maptools)
data(wrld_simpl)
plot(wrld_simpl, add = TRUE)

ファイルからディメンション情報を読み取ることでこれを実現する方法は他にもたくさんありますが、これは少なくとも、ハードワークのほとんどを実行するためのショートカットです。

ラスターからの画像のプロット

于 2012-12-04T10:35:59.900 に答える
7

@mdsumnerのはるかに優れたソリューションを補完するのと同じように、ライブラリncdfのみを使用してこれを行う方法があります。

library(ncdf)
f <- open.ncdf("LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040101.nc")
A <- get.var.ncdf(nf,"soil_moisture_c")

必要なのは、コヒーレントなx軸とy軸を持つために寸法を見つけることだけです。NetCDFオブジェクトのディメンションを見ると、次のように表示されます。

str(f$dim)
List of 2
 $ Longitude:List of 8
  ..$ name         : chr "Longitude"
  ..$ len          : int 1440
  ..$ unlim        : logi FALSE
  ..$ id           : int 1
  ..$ dimvarid     : num 2
  ..$ units        : chr "degrees_east"
  ..$ vals         : num [1:1440(1d)] -180 -180 -179 -179 -179 ...
  ..$ create_dimvar: logi TRUE
  ..- attr(*, "class")= chr "dim.ncdf"
 $ Latitude :List of 8
  ..$ name         : chr "Latitude"
  ..$ len          : int 720
  ..$ unlim        : logi FALSE
  ..$ id           : int 2
  ..$ dimvarid     : num 1
  ..$ units        : chr "degrees_north"
  ..$ vals         : num [1:720(1d)] 89.9 89.6 89.4 89.1 88.9 ...
  ..$ create_dimvar: logi TRUE
  ..- attr(*, "class")= chr "dim.ncdf"

したがって、寸法は次のとおりです。

 f$dim$Longitude$vals -> Longitude
 f$dim$Latitude$vals -> Latitude

Latitudeこれで、逆ではなく90から-90になります。これにより、次のimageようになります。

 Latitude <- rev(Latitude)
 A <- A[nrow(A):1,]

最後に、お気づきのとおり、オブジェクトAのxとyが反転しているため、転置する必要があります。NA値は、何らかの理由で値で表されます-32767。

A[A==-32767] <- NA
A <- t(A)

そして最後にプロット：

image(Longitude, Latitude, A)
library(maptools)
data(wrld_simpl)
plot(wrld_simpl, add = TRUE)

ここに画像の説明を入力してください

編集： 31個のファイルでこれを行うには、ファイル名のベクトルとそれらを保存したディレクトリを呼び出しますncfiles（yourpath簡単にするために、変数は常に呼び出されsoil_moisture_c、NAは常に呼び出されると仮定します-32767）：

ncfiles
 [1] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040101.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040102.nc"
 [3] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040103.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040104.nc"
 [5] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040105.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040106.nc"
 [7] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040107.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040108.nc"
 [9] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040109.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040110.nc"
[11] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040111.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040112.nc"
[13] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040113.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040114.nc"
[15] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040115.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040116.nc"
[17] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040117.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040118.nc"
[19] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040119.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040120.nc"
[21] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040121.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040122.nc"
[23] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040123.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040124.nc"
[25] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040125.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040126.nc"
[27] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040127.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040128.nc"
[29] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040129.nc" "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040130.nc"
[31] "LPRM-AMSR_E_L3_D_SOILM3_V002-20120520T173800Z_20040131.nc"

yourpath
 [1] "C:\\Users"

library(ncdf)
library(maptools)
data(wrld_simpl)
for(i in 1:length(ncfiles)){
    f <- open.ncdf(paste(yourpath,ncfiles[i], sep="\\"))
    A <- get.var.ncdf(f,"soil_moisture_c")
    f$dim$Longitude$vals -> Longitude
    f$dim$Latitude$vals -> Latitude
    Latitude <- rev(Latitude)
    A <- A[nrow(A):1,]
    A[A==-32767] <- NA
    A <- t(A)
    close.ncdf(f) # this is the important part
    png(paste(gsub("\\.nc","",ncfiles[i]), "\\.png", sep="")) # or any other device such as pdf, jpg...
    image(Longitude, Latitude, A)
    plot(wrld_simpl, add = TRUE)
    dev.off()
    }
于 2012-12-04T10:59:10.610 に答える
1

CDOを使用して、コマンドラインから緯度を単純に反転することもできます。

cdo invertlat file.nc file_inverted.nc
于 2017-09-01T15:35:07.057 に答える