Project

General

Profile

Download (3.44 KB) Statistics
| Branch: | Revision:
1
# R script for creating a plot of mean elevation profile containing four lines 
2
# (aster2, srtm, uncorrected fused, blend gaussian fused) for each nine region
3
# The results (plots) are saved at 
4
#   "/data/project/organisms/DEM?Yuni/documents/meanElv"
5
#
6
# Original author: Jim Regetz
7
# [13-JUL-2011]
8
#   Edits by Yuina to auto-apply for nine regions.
9
#   
10
#   9-Dec-2011
11
#   Yuina Nunokawa 
12

    
13

    
14

    
15
library(raster)
16

    
17
aster.dir      <- "/data/project/organisms/DEM/Yuni/Data/aster2"
18
srtm.dir       <- "/data/project/organisms/DEM/Yuni/Data/srtm"
19
aster.straddle <- list.files(aster.dir, pattern="^aster2.*_straddle.tif")
20
srtm.below     <- list.files(srtm.dir, pattern="^srtm_.*_below.tif") 
21
fused.straddle <- list.files(aster.dir, pattern="^fused.*_straddle.tif")
22
fused.bgau     <- list.files(aster.dir, pattern="^fused.*_blendgau.tif")
23

    
24

    
25

    
26
mElevation <- function(d.aster2, d.srtm, d.uncor2, d.bg2, file.name, ylim, aster.dir, srtm.dir){
27

    
28
    # extract raster latitudes for later
29
    lats2400 <- yFromRow(d.aster2, 1:nrow(d.aster2))
30
    lats1200 <- yFromRow(d.srtm, 1:nrow(d.srtm))
31

    
32
    # initialize output pdf device driver
33
    name <- paste("mElevation_", file.name, ".pdf", sep="") 
34
    pdf(file.path("/data/project/organisms/DEM/Yuni/documents/meanElv", file=name), height=8, width=11.5)
35
    
36
    par(mfrow=c(2,2), omi=c(1,1,1,1))
37
    ylim <- ylim
38
   
39
    # plot latitudinal profiles of mean elevation for aster2 and srtm 
40
    plot(lats2400, rowMeans(as.matrix(d.aster2), na.rm=TRUE), type="l",
41
        xlab="Latitude", ylab="Mean elevation", ylim=ylim)
42
    text(min(lats1200), min(ylim)+0.5, pos=4, font=3, labels="AGDEM2")
43
    lines(lats1200, rowMeans(as.matrix(d.srtm), na.rm=TRUE), col="blue")
44
    legend("bottomright", legend=c("SRTM", "ASTER2"), col=c("blue",
45
        "black"), lty=c(1, 1), bty="n")
46
    abline(v=60, col="red", lty=2)
47
    mtext(expression(paste("Latitudinal profiles of mean elevation", file.name)), adj=0, line=2, font=2)
48
    
49

    
50
    # plot latitudial profiles of mean elevation for Uncorrected fused and blend gaussian fused
51
    plot(lats2400, rowMeans(as.matrix(d.uncor2), na.rm=TRUE), type="l",
52
        xlab="Latitude", ylab="Mean elevation", ylim=ylim)
53
    text(min(lats2400), min(ylim)+0.5, pos=4, font=3, labels="simple fuse and gaussian")
54
    lines(lats2400, rowMeans(as.matrix(d.bg2), na.rm=TRUE), col="red")
55
    legend("bottomright", legend=c("Simple Fused", "Gaussian Blend"), col=c("black","red"), lty=c(1, 1), bty="n")
56
    abline(v=60, col="red", lty=2)
57
    dev.off () 
58

    
59
}
60

    
61

    
62

    
63

    
64

    
65
num.regions <- 9
66
plotElevation <-lapply(1:num.regions, function(region, aster.dir, srtm.dir){
67
    d.aster2 <- raster(file.path(aster.dir, aster.straddle[region])) 
68
    d.srtm <- raster(file.path(srtm.dir, srtm.below[region]))  
69
    d.uncor2 <- raster(file.path(aster.dir, fused.straddle[region]))
70
    d.bg2 <- raster(file.path(aster.dir, fused.bgau[region])) 
71
    file.name <- substr(srtm.below[region], 6, 13)
72
    
73
    if (region==1){
74
        ylim <- c(120, 150)
75
    }else if(region==2){
76
        ylim <- c(140, 175)
77
    }else if (region==3){   
78
        ylim <- c(350, 550)
79
    }else if (region==4){
80
        ylim <-c(50, 450)
81
    }else if (region==5){
82
        ylim <- c(60, 200)
83
    }else if (region==6){
84
        ylim <- c(0, 90)
85
    }else if (region==7){
86
        ylim <- c(60, 100)
87
    }else if (region==8){
88
        ylim <- c(600, 750)
89
    }else{
90
        ylim <- c(30, 280)
91
    }
92

    
93
    mElevation(d.aster2, d.srtm, d.uncor2, d.bg2, file.name, ylim, aster.dir, srtm.dir)
94
}, aster.dir=aster.dir, srtm.dir=srtm.dir)
95

    
(14-14/23)