计算并绘制任意 rasterLayer 的矢量场

2023-12-23

问题陈述:

With ggquiver::geom_quiver()只要我们知道,我们就可以绘制向量场x, y, xend, and yend.

  1. 我如何计算任意的这些参数RasterLayer海拔?
  2. 如何确保这些箭头的大小指示该特定向量的斜率,以便箭头显示与该位置处的梯度成比例的不同长度(例如,下面的第一张图)?

背景:

# ggquiver example

library(tidyverse)
library(ggquiver)
expand.grid(x=seq(0,pi,pi/12), y=seq(0,pi,pi/12)) %>%
  ggplot(aes(x=x,y=y,u=cos(x),v=sin(y))) +
  geom_quiver()

相关方法使用rasterVis::vectorplot,这依赖于raster::terrain(假设场单位 == CRS 单位)计算并绘制矢量场。源代码在这里 https://github.com/oscarperpinan/rastervis/blob/master/R/vectorplot.R.

library(raster)
library(rasterVis)
r <- getData('alt', country='FRA', mask=TRUE)
r <- aggregate(r, 20)
vectorplot(r, par.settings=RdBuTheme())

结论:

为了审查,我想采取任意rasterLayer的高程,将其转换为data.frame,计算x, y, xmax, and ymax高程矢量场的分量,用于调整箭头的大小,以便它们显示该点的相对坡度(如上面的图 1 和图 2 所示),并用ggquiver。就像是:

names(r) <- "z"
rd <- as.data.frame(r, xy=TRUE)

# calculate x, y, xend, yend for gradient vectors, add to rd, then plot

ggplot(rd) + 
  geom_raster(aes(x, y, fill = z)) + 
  geom_quiver(aes(x, y, xend, yend))

实际上,您要求的是转换 2D标量场 https://en.wikipedia.org/wiki/Scalar_field into a 矢量场 https://en.wikipedia.org/wiki/Vector_field。有几种不同的方法可以做到这一点。

光栅包包含函数terrain,这会创建新的栅格图层,该图层将为您提供每个点处所需矢量的角度(即aspect)及其大小(slope)。我们可以使用一点三角学将它们转换为南北和东西向基本向量ggquiver并将它们添加到我们的原始栅格中,然后将整个数据框转换为数据框。*

terrain_raster <- terrain(r, opt = c('slope', 'aspect'))
r$u <- terrain_raster$slope[] * sin(terrain_raster$aspect[])
r$v <- terr$slope[] * cos(terr$aspect[])
rd <- as.data.frame(r, xy = TRUE)

然而,在大多数情况下,这不会产生好的情节。如果您不首先聚合栅格,则图像上的每个像素都会有一个渐变,这将无法很好地绘制。另一方面,如果你do聚合后,您将拥有一个不错的矢量场,但您的栅格将看起来“块状”。因此,为绘图使用单个数据框可能不是最好的方法。

以下函数将获取一个栅格并用叠加的矢量场绘制它。您可以在不影响栅格的情况下调整矢量场的聚合程度,并且可以为栅格指定任意颜色矢量。

raster2quiver <- function(rast, aggregate = 50, colours = terrain.colors(6))
{
  names(rast) <- "z"
  quiv <- aggregate(rast, aggregate)
  terr <- terrain(quiv, opt = c('slope', 'aspect'))
  quiv$u <- terr$slope[] * sin(terr$aspect[])
  quiv$v <- terr$slope[] * cos(terr$aspect[])
  quiv_df <- as.data.frame(quiv, xy = TRUE)
  rast_df <- as.data.frame(rast, xy = TRUE)

  print(ggplot(mapping = aes(x = x, y = y, fill = z)) + 
          geom_raster(data = rast_df, na.rm = TRUE) + 
          geom_quiver(data = quiv_df, aes(u = u, v = v), vecsize = 1.5) +
          scale_fill_gradientn(colours = colours, na.value = "transparent") +
          theme_bw())

  return(quiv_df)
}

因此,在您的法国示例中尝试一下,在首先定义类似的调色板之后,我们得到

pal <- c("#B2182B", "#E68469", "#D9E9F1", "#ACD2E5", "#539DC8", "#3C8ABE", "#2E78B5")

raster2quiver(getData('alt', country = 'FRA', mask = TRUE), colours = pal)

现在要证明它可以在随意的栅格(假设它已指定投影)让我们在该图像上测试它(转换为栅格)。这次,我们的分辨率较低,因此我们选择较小的合计值。我们还将为最低值选择透明颜色,以提供更好的绘图:

rast <- raster::raster("https://i.stack.imgur.com/tXUXO.png")

# Add a fake arbitrary projection otherwise "terrain()" doesn't work:
projection(rast) <- "+proj=lcc +lat_1=48 +lat_2=33 +lon_0=-100 +ellps=WGS84"

raster2quiver(rast, aggregate = 20, colours = c("#FFFFFF00", "red"))

* I should point out that geom_quiver's mapping aesthetic takes arguments called u and v, which represent the basis vectors pointing North and East. The ggquiver package converts them to xend and yend values using stat_quiver. If you prefer to use xend and yend values you could just use geom_segment to plot your vector field, but this makes it more complex to control the appearance of the arrows. Hence, this solution will find the magnitude of the u and v values instead.

本文内容由网友自发贡献,版权归原作者所有,本站不承担相应法律责任。如您发现有涉嫌抄袭侵权的内容,请联系:hwhale#tublm.com(使用前将#替换为@)

计算并绘制任意 rasterLayer 的矢量场 的相关文章

  • 将绘图调用拆分为多个块

    我正在编写一个图的解释 其中我基本上将在第一个块中创建图 然后描述该输出 并在第二个块中添加一个轴 然而 似乎每个块都会强制一个新的绘图环境 因此当我们尝试使用以下命令运行块时会出现错误axis独自的 观察 output html docu
  • 在 R 中绘制 Likert 变量的堆积条形图

    假设我有一个如下所示的数据框 P Q1 Q2 1 1 4 1 2 2 3 4 3 1 1 4 其中的列告诉我哪个人相应地回答了问题 q1 q2 中的哪一个 这些问题需要按照 4 分李克特量表进行回答 例如 批准 表示 1 稍微批准 表示 2
  • twitterR 和 ROAuth R 软件包安装

    我在安装 CRAN 上的 twitteR 和 RAOuth 软件包时遇到一些问题 我尝试了几种不同的方法 在 Windows 下使用源代码 在 Ubuntu 下使用 RStudio 我尝试了以下命令 sudo apt get install
  • 朴素贝叶斯分类器仅基于先验概率做出决策

    我试图根据推文的情绪将推文分为三类 买入 持有 卖出 我正在使用 R 和包 e1071 我有两个数据框 一个训练集和一组需要预测情绪的新推文 训练集数据框 text sentiment this stock is a good buy Bu
  • 只读取选定的列

    谁能告诉我如何仅读取下面每年数据的前 6 个月 7 列 例如使用read table Year Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec 2009 41 27 25 31 31 39 2
  • 相当于 min() 的 rowMeans()

    我在 R 邮件列表上多次看到这个问题 但仍然找不到满意的答案 假设我有一个矩阵m m lt matrix rnorm 10000000 ncol 10 我可以通过以下方式获得每行的平均值 system time rowMeans m use
  • 文本挖掘 pdf 文件/词频问题

    我正在尝试挖掘一篇具有丰富 pdf 编码和图表的文章的 pdf 我注意到 当我挖掘一些 pdf 文档时 我得到的高频词是 phi taeoe toe sigma gamma 等 它与某些 pdf 文档配合良好 但与其他文档配合使用时却得到这
  • 旋转 Markdown 的表格 pdf 输出

    我想将 pdf 上的表格输出旋转 90 度 我正在使用 Markdown 生成报告并kable循环显示表格 如果可以的话我想继续使用kable因为还有很多其他依赖于它的东西我没有包含在这个 MWE 中 这是一个简单的例子 使用iris数据集
  • 如何仅删除单括号并保留配对的括号

    你好 我亲爱的老师 R 用户朋友们 我最近开始认真学习正则表达式 最近我遇到了一种情况 我们只想保留配对括号 并省略未配对的 这是我的样本数据 structure list t1 c Book Pg 1 Website Online Jou
  • 如何从 R keras 中的类似生成器的数据中评估()和预测()

    我有以下代码 数据集可以下载here https www dropbox com s qjt5o31oyqj10m8 data tar gz dl 0 or here https www kaggle com c dogs vs cats
  • 计算 R 中各列的唯一值

    我正在尝试创建一个新变量 其中包含来自两个不同列的字符串值的唯一计数 所以我有这样的东西 例如 A tibble 4 x 2 names partners
  • 在 Rcpp 中使用其他包中的 C 函数

    我试图从 C 函数中的 cubature 包调用 C 例程来执行多维积分 我试图重现的基本 R 示例是 library cubature integrand lt function x sin x adaptIntegrate integr
  • 所有 x 轴标签未以 45 度显示

    I m having the code as like below But I m not getting all the x axis labels and it is not displaying in 45 degree when I
  • 为什么 R 更新后 sim_slopes() 中会出现此错误?

    我正在尝试使用 交互 包来创建简单斜率的约翰逊 尼曼图 但是 当尝试运行 sim slopes 函数时 出现以下错误 直到我将R更新到4 2 2 我才没有遇到这个问题 我使用的是 macOS Ventura 13 1 Error class
  • python 相当于 R 中的 get() (= 使用字符串检索符号的值)

    在 R 中 get s 函数检索名称存储在字符变量 向量 中的符号的值s e g X lt 10 r lt XVI s lt substr r 1 1 X get s 10 取罗马数字的第一个符号r并将其转换为其等效整数 尽管花了一些时间翻
  • SPSS 中的标准化残差与 R rstandard(lm()) 不匹配

    在寻找 R 相关解决方案时 我发现 R 和 SPSS 版本 24 在计算简单线性模型中的标准化残差方面存在一些不一致 看来SPSS所谓的标准化残差匹配 R学生化残差 我完全不认为某处存在软件错误 但显然这两个程序之间存在差异 看看这个例子
  • Shiny:动态数据框构建; renderUI、观察、reactiveValues

    我认为如何使用 Shiny 的 renderUI 功能动态子集数据的问题经常出现 但我很难理解何时使用 renderUI 带有 uiOutput 而不是其他功能 包括观察 反应 反应值甚至条件面板 我想构建一个完全交互式的数据框架 其中每个
  • 麦当劳 omega:R 中的警告

    我正在计算几种不同尺度的欧米茄 并在 R 中使用不同的 omega 函数获取不同比例的不同警告消息 我的问题是如何解释这些警告以及报告检索到的 omega 统计数据是否安全 当我使用 从 alpha 到 omega 内部一致性估计普遍问题的
  • 如何将plot中的单变量列表图表转换为ggplot2格式?

    我正在搜索 但仍然找不到一个非常简单的问题的答案 我们如何使用 R 中的 ggplot2 生成一个变量的简单线图 我正在分析时间序列数据 并且想要对图表进行更复杂的操作 我认为如果我使用 ggplot2 代替会更好plot It works
  • 在 Shiny 中的用户会话之间共享反应数据集

    我有一个相当大的反应数据集 该数据集是通过轮询文件然后按预定义的时间间隔读取该文件而派生的 数据更新频繁 需要不断重新加载 诚然 重新加载可以增量完成并附加到 R 中的现有对象 但事实并非如此 然而目前 尽管会话中的数据相同 但此操作是针对

随机推荐

  • Java 中接口的重要性[关闭]

    很难说出这里问的是什么 这个问题是含糊的 模糊的 不完整的 过于宽泛的或修辞性的 无法以目前的形式得到合理的回答 如需帮助澄清此问题以便重新打开 访问帮助中心 help reopen questions 假设我们有两个班级 Tiger an
  • Windows Remote 上 OpenGL 的现状和解决方案 [关闭]

    Closed 这个问题不符合堆栈溢出指南 help closed questions 目前不接受答案 OpenGL 和 Windows Remote 不能很好地配合 此问题的解决方案取决于用例 并且答案分散在网络的各个角落 当我开始研究这个
  • 在 Crystal Reports 中将行数据转置为列

    我从存储过程返回以下数据 Staff Category Amount Bob Art 123 Bob Sport 777 Bob Music 342 Jeff Art 0 Jeff Sport 11 Jeff Music 27 即使金额为零
  • Mac OS X:从目录服务获取当前用户的当前用户名和主目录

    我的应用程序是用 Objective C 编写的 如何通过目录服务获取当前登录用户的用户名和主目录 细节 我的 Cocoa 应用程序使用 getenv USER getenv HOME 获取当前用户名和主目录 显然 如果用户通过目录服务登录
  • leafletjs 添加可滚动弹出窗口?

    使用带有弹出窗口的 leafletjs 当我的弹出窗口包含最少的文本时 一切正常 如果我把它们做得更大 它们仍然可以正常工作 如果我添加太多 我会添加 maxHeight 到弹出窗口 它使弹出窗口可滚动 如果我一开始没有足够的内容来填充页面
  • 如何在 MediaPlayer setDataSource 中包含 http 标头?

    我正在将 URI 传递给设置数据源方法 http developer android com reference android media MediaPlayer html setDataSource 28java lang String
  • 使用动态发出的 POCO 进行快速序列化和反序列化

    我目前正在将 SQL 表行序列化为二进制格式 以实现高效存储 我将二进制数据序列化 反序列化为List每行 我正在尝试升级它以使用 POCO 它将动态生成 发出 每列一个字段 我在网上搜索了几个小时 偶然发现了像 EF T4 Expando
  • 如何替换字符串中的最后一个单词

    有谁知道如何替换字符串中的最后一个单词 目前我正在做 someStr someStr replace someStr substring someStr lastIndexOf 1 New Word 上面的代码替换了字符串中出现的每个单词
  • 如何使用 Flexbox 实现浮动侧边栏布局,内容环绕侧边栏?

    我正在尝试实现在桌面上看起来像这样的响应式布局 在手机上就像这样 请注意以下要求 侧边栏应仅占据适合内容所需的垂直空间 在侧边栏下方 主要部分的内容应占据整个宽度 在移动设备上 侧边栏应显示在主要内容下方 这是一个包含我最初的 HTML 和
  • 使用装饰器自动注册类方法

    我希望能够创建一个 python 装饰器 自动在全局存储库中 注册 类方法 带有一些属性 示例代码 class my class object register prop1 prop2 def my method arg1 arg2 met
  • HTML5 视频 MEDIA_ERR_DECODE 随机发生

    我正在开发一个包含 6 个音频和视频元素的项目 这些元素依次播放 发出前的代码顺序是这样的 预加载所有媒体资源直到 canplaythrough 播放视频 1 停止 video 1 并播放 audio 1 停止音频 1 并再次播放视频 1
  • 使用 sed/awk 替换文本文件的部分

    我正在尝试替换文件中 begin 和 end 之间的文本 如下所示 begin block param1 param2 end 读作 begin block param1 value1 param2 value2 end 基本上 我取消注释
  • 从 JDialog 返回值; dispose()、setVisible(false) - 示例

    我知道 这个问题在 SO 中经常出现 比如here https stackoverflow com questions 4089311 how can i return a value from a jdialog box to the p
  • 如何在 rechart 中将工具提示放置在条形图的顶部?

    Problem 我创建了一个带有自定义工具提示的条形图 现在我需要的是将工具提示放置在栏的顶部而不是图表区域内 就像这张照片一样 这就是现在的样子 在这里我提供了我如何组织我的代码 import React Component from r
  • 如何将一个计算列的每一行除以另一个计算列的总和?

    我无法使用此示例数据得到正确的除法 Calculated column Another calc column 48 207 257 370 518 138 489 354 837 478 1 005 648 1 021 2 060 1 4
  • 如何获得 COUNT(column) ... GROUP BY 来使用索引?

    我有一个表 col1 col2 其索引为 col1 col2 该表中有数百万行 我想运行一个查询 SELECT col1 COUNT col2 WHERE col1 NOT IN
  • 在 Highcharts/Highmaps 中的数据集之间切换

    我正在尝试按县创建美国生猪数量的分区统计图 我有六个不同的数据集 1987 1992 1997 2002 2007 2012 代表每五年进行一次的美国农业部农业普查 截至目前 每个数据集都位于其自己的 Highmaps 索引文件中 2012
  • 更漂亮:在函数和括号之间添加空格

    在 VSCode 中 每次保存 JS 文件时 Prettier删除 function 关键字及其括号之间的空格 它改变了这一点 function parameter To this function parameter 但我想保留空间 有些
  • 停止设备的当前密码要求

    我的网站的编辑用户部分遇到问题 由于某种原因 我在尝试编辑用户时不断收到错误 当前密码不能为空 我们使用 devise 来管理用户 但我似乎无法在任何地方找到会生成此错误的代码 这是表单的代码 semantic form for resou
  • 计算并绘制任意 rasterLayer 的矢量场

    问题陈述 With ggquiver geom quiver 只要我们知道 我们就可以绘制向量场x y xend and yend 我如何计算任意的这些参数RasterLayer海拔 如何确保这些箭头的大小指示该特定向量的斜率 以便箭头显示