2014-04-01 4 views
-1

여기에서 다운로드 할 수있는 글로벌 면적 나타내는 209091 요소 1D 이진 배열이 포함 된 파일이 있습니다 내가 제공하는 보조 행과 열 파일을 사용하여 1 차원 데이터 배열에서 풀을 만들려면 ftp://sidads.colorado.edu/DATASETS/nsidc0451_AMSRE_Land_Parms_v01/AMSRE_flags_2002/ 을 .globland_r ftp://sidads.colorado.edu/DATASETS/nsidc0451_AMSRE_Land_Parms_v01/AMSRE_ancil/R을 사용하여 1D 배열에서 그리드를 만드는 방법?

는이 목적을 위해 matlab에 작성된 코드가 나는 R이 MATLAB 코드를 번역 할하지만 다음 matlab에에게

 function [gridout, EASE_r, EASE_s] = mkgrid_global(x) 
    %MKGRID_GLOBAL(x) Creates a matrix for mapping 
    % gridout = mkgrid_global(x) uses the 2090887 element array (x) and returns 

    %Load ancillary EASE grid row and column data, where <MyDir> is the path to 
    %wherever the globland_r and globland_c files are located on your machine. 
    fid = fopen('C:\MyDir\globland_r','r'); 
    EASE_r = fread(fid, 209091, 'int16'); 
    fclose(fid); 


    fid = fopen('C:\MyDir\globland_c','r'); 
    EASE_s = fread(fid, 209091, 'int16'); 
    fclose(fid); 



    gridout = NaN.*zeros(586,1383); 

    %Loop through the elment array 
    for i=1:1:209091  

    %Distribute each element to the appropriate location in the output 
    %matrix (but MATLAB is 
    %(1,1) 

    end 

편집을 모르는 : 여기에서 다운로드 할 수있는 globland_c 티 @mdsumner의 그 해결책 :

파일 MLLATLSBMLLONLSB (4-byte integers)을 위해 위도와 경도 (multiply by 1e-5)를 포함하는 지리적 위치 전체 글로벌 EASE grid matrix (586×1383) MLLATLSBMLLONLSB 여기에서 다운로드 할 수 있습니다 ftp://sidads.colorado.edu/DATASETS/nsidc0451_AMSRE_Land_Parms_v01/AMSRE_ancil/

답변

0
## the sparse dims, literally the xcol * yrow indexes 
dims <- c(1383, 586) 

cfile <- "ftp://sidads.colorado.edu/DATASETS/nsidc0451_AMSRE_Land_Parms_v01/AMSRE_ancil/globland_c" 
rfile <- "ftp://sidads.colorado.edu/DATASETS/nsidc0451_AMSRE_Land_Parms_v01/AMSRE_ancil/globland_r" 

## be nice, don't abuse this 
col <- readBin(cfile, "integer", n = prod(dims), size = 2, signed = FALSE) 
row <- readBin(rfile, "integer", n = prod(dims), size = 2, signed = FALSE) 

## example data file 
fdat <- "ftp://sidads.colorado.edu/DATASETS/nsidc0451_AMSRE_Land_Parms_v01/AMSRE_flags_2002/flags_2002170A.bin" 

dat <- readBin(fdat, "integer", n = prod(dims), size = 1, signed = FALSE) 


## now get serious 
m <- matrix(as.integer(NA), dims[2L], dims[1L]) 
m[cbind(row + 1L, col + 1L)] <- dat 


image(t(m)[,dims[2]:1], col = rainbow(length(unique(m)), alpha = 0.5)) 

어쩌면 우리가 할 수있는 이지도 투영도 재구성하십시오.

flon <- "MLLONLSB" 
flat <- "MLLATLSB" 
## the key is that these are integers, floats scaled by 1e5 
lon <- readBin(flon, "integer", n = prod(dims), size = 4) * 1e-5 
lat <- readBin(flat, "integer", n = prod(dims), size = 4) * 1e-5 

## this is all we really need from now on 
range(lon) 
range(lat) 

library(raster) 
library(rgdal) ## need for coordinate transformation 
ex <- extent(projectExtent(raster(extent(range(lon), range(lat)), crs = "+proj=longlat"), "+proj=cea")) 

grd <- raster(ncols = dims[1L], nrows = dims[2L], xmn = xmin(ex), xmx = xmax(ex), ymn = ymin(ex), ymx = ymax(ex), crs = "+proj=cea") 

거기에는 "out half pixel"오류가있을 수 있습니다.

테스트

plot(setValues(grd, m), col = rainbow(max(m, na.rm = TRUE), alpha = 0.5)) 

Hohum

library(maptools) 
data(wrld_simpl) 
plot(spTransform(wrld_simpl, CRS(projection(grd))), add = TRUE) 
우리는 지금 다음, 우리의 "GRD"템플릿과 일치하는 특정 DAT 파일을 읽고 단지 그 값으로 템플릿을 채우기 위해 유효한 cellnumbers을 절약 할 수 있습니다

cellnumbers에 따라. 더구나, 누군가는이 경로를 일찌기 거의 밟았지만 많이 얻은 것 같지 않습니다.

How to identify lat and long for a global matrix?

관련 문제