Fortran图形绘制与可视化

Fortran 图形绘制基础

在 Fortran 中进行图形绘制并非其原生强项，但借助一些第三方库，我们能够实现强大的图形绘制与可视化功能。其中，较为常用的是 PGPLOT 库。PGPLOT 是一个便携式的 Fortran 绘图库，可用于生成各种类型的二维和三维图形。

安装 PGPLOT 库

1. 下载源文件：从 PGPLOT 的官方网站（http://www.astro.caltech.edu/~tjp/pgplot/）下载适合你系统的源文件压缩包。例如，对于类 Unix 系统，可能会下载 pgplot5.2.tar.gz。
2. 解压与编译：解压下载的压缩包，进入解压后的目录，通常会看到一系列的源文件和 Makefile。根据你的系统类型，可能需要对 Makefile 进行一些调整。例如，在 Linux 系统下，如果你的编译器是 GCC，可能需要确保 F77 和 CC 变量设置为正确的 GCC 版本路径。然后执行 make 命令进行编译。编译成功后，会生成一些库文件（如 libpgplot.a）和可执行文件。
3. 安装：将编译生成的库文件和头文件安装到系统的相关目录中，以便 Fortran 编译器能够找到它们。例如，你可以将库文件拷贝到 /usr/local/lib，头文件拷贝到 /usr/local/include，并确保这些目录在编译器的搜索路径中。

简单二维图形绘制

下面通过一个简单的示例展示如何使用 PGPLOT 库在 Fortran 中绘制二维线图。

program plot_example
    use pgplot
    implicit none
    integer :: i, ndata
    real :: x(100), y(100)
    ndata = 100
   ! 生成数据
    do i = 1, ndata
        x(i) = real(i)
        y(i) = sin(x(i) / 10.0)
    end do
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('/XWINDOW')
   ! 设置绘图区域
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv(0.0, real(ndata), -1.2, 1.2, 0, 0)
   ! 绘制线图
    call pgline(ndata, x, y)
   ! 添加标题和坐标轴标签
    call pglab('X 轴', 'Y 轴', '正弦函数曲线')
   ! 关闭 PGPLOT
    call pgclos()
end program plot_example

1. 数据生成：在上述代码中，首先定义了数组 x 和 y 来存储数据点。通过循环生成了 x 从 1 到 100 的值，并计算对应的 y = sin(x / 10.0) 值。
2. PGPLOT 初始化：call pgpap(0.0, 0.0) 用于设置画笔参数，这里使用默认值。call pgopen('/XWINDOW') 打开绘图设备，这里选择了 X Window 作为输出设备。如果在不同的系统或希望输出到文件（如 PostScript 文件），可以更改设备名称。
3. 绘图区域设置：call pgsvp(0.15, 0.95, 0.15, 0.95) 设置了绘图区域在整个图形窗口中的比例，即从窗口的 15% 位置开始到 95% 位置结束。call pgenv(0.0, real(ndata), -1.2, 1.2, 0, 0) 设置了坐标轴的范围，x 轴从 0 到 ndata，y 轴从 -1.2 到 1.2。
4. 图形绘制：call pgline(ndata, x, y) 用于绘制线图，将 x 和 y 数组中的数据点用线段连接起来。
5. 添加标签：call pglab('X 轴', 'Y 轴', '正弦函数曲线') 为图形添加了坐标轴标签和标题。
6. 关闭设备：最后，call pgclos() 关闭绘图设备，释放资源。

高级二维图形绘制

绘制散点图

散点图在展示数据分布关系方面非常有用。以下是使用 PGPLOT 绘制散点图的示例代码。

program scatter_plot
    use pgplot
    implicit none
    integer :: i, ndata
    real :: x(200), y(200)
    ndata = 200
   ! 生成随机数据
    call random_seed()
    do i = 1, ndata
        call random_number(x(i))
        call random_number(y(i))
    end do
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('/XWINDOW')
   ! 设置绘图区域
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv(0.0, 1.0, 0.0, 1.0, 0, 0)
   ! 绘制散点图
    call pgsls(1)
    call pgsci(2)
    call pgpt(ndata, x, y, 1)
   ! 添加标题和坐标轴标签
    call pglab('X 轴', 'Y 轴', '随机散点图')
   ! 关闭 PGPLOT
    call pgclos()
end program scatter_plot

1. 数据生成：通过 call random_seed() 初始化随机数种子，然后在循环中使用 call random_number 生成 0 到 1 之间的随机数填充 x 和 y 数组。
2. 绘图设置：与之前类似，先初始化 PGPLOT，设置绘图区域。这里通过 call pgsls(1) 设置线条样式为实线（如果绘制散点图的连线），call pgsci(2) 设置颜色索引为 2（具体颜色取决于设备，通常为红色）。call pgpt(ndata, x, y, 1) 用于绘制散点，其中 1 表示使用实心圆作为点的符号。

绘制柱状图

柱状图常用于比较不同类别数据的大小。以下是绘制柱状图的代码示例。

program bar_chart
    use pgplot
    implicit none
    integer :: i, ncategories
    real :: x(5), height(5)
    ncategories = 5
   ! 定义类别和高度数据
    x = [1.0, 2.0, 3.0, 4.0, 5.0]
    height = [12.0, 18.0, 9.0, 20.0, 15.0]
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('/XWINDOW')
   ! 设置绘图区域
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv(0.5, 5.5, 0.0, 25.0, 0, 0)
   ! 绘制柱状图
    do i = 1, ncategories
        call pgsls(1)
        call pgsci(3)
        call pgbars(x(i), 1.0, height(i), 0.0, 1)
    end do
   ! 添加标题和坐标轴标签
    call pglab('类别', '数值', '柱状图')
   ! 关闭 PGPLOT
    call pgclos()
end program bar_chart

1. 数据定义：数组 x 定义了类别位置，height 数组定义了每个类别的高度。
2. 绘图操作：在循环中，使用 call pgbars(x(i), 1.0, height(i), 0.0, 1) 绘制柱状图。其中 x(i) 是柱子的中心位置，1.0 是柱子的宽度，height(i) 是柱子的高度，0.0 是柱子底部的偏移（这里为 0），1 表示绘制填充的柱子。

三维图形绘制

使用 PGPLOT 绘制三维表面图

PGPLOT 也支持三维图形绘制，下面是绘制三维表面图的示例。

program surface_plot
    use pgplot
    implicit none
    integer :: i, j, nx, ny
    real :: x(100), y(100), z(100, 100)
    nx = 100
    ny = 100
   ! 生成数据
    do i = 1, nx
        x(i) = (i - 1) / real(nx - 1) * 2.0 - 1.0
    end do
    do j = 1, ny
        y(j) = (j - 1) / real(ny - 1) * 2.0 - 1.0
    end do
    do i = 1, nx
        do j = 1, ny
            z(i, j) = exp(-(x(i) ** 2 + y(j) ** 2))
        end do
    end do
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('/XWINDOW')
   ! 设置三维绘图
    call pgslw(2)
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv3d(-1.0, 1.0, -1.0, 1.0, 0.0, 1.0, 1)
    call pgsurf(nx, ny, x, y, z)
   ! 添加标题和坐标轴标签
    call pglab3d('X 轴', 'Y 轴', 'Z 轴', '三维表面图')
   ! 关闭 PGPLOT
    call pgclos()
end program surface_plot

1. 数据生成：通过嵌套循环生成 x 和 y 数组的值，范围从 -1 到 1。然后根据 x 和 y 计算 z 数组的值，这里使用了高斯函数 z = exp(-(x ** 2 + y ** 2))。
2. 三维绘图设置：call pgslw(2) 设置线条宽度为 2。call pgenv3d(-1.0, 1.0, -1.0, 1.0, 0.0, 1.0, 1) 设置三维坐标轴的范围。call pgsurf(nx, ny, x, y, z) 用于绘制三维表面图。

绘制三维散点图

program scatter_3d_plot
    use pgplot
    implicit none
    integer :: i, ndata
    real :: x(300), y(300), z(300)
    ndata = 300
   ! 生成随机数据
    call random_seed()
    do i = 1, ndata
        call random_number(x(i))
        call random_number(y(i))
        call random_number(z(i))
    end do
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('/XWINDOW')
   ! 设置三维绘图区域
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv3d(0.0, 1.0, 0.0, 1.0, 0.0, 1.0, 1)
   ! 绘制三维散点图
    call pgsls(1)
    call pgsci(4)
    call pgpt3d(ndata, x, y, z, 1)
   ! 添加标题和坐标轴标签
    call pglab3d('X 轴', 'Y 轴', 'Z 轴', '三维随机散点图')
   ! 关闭 PGPLOT
    call pgclos()
end program scatter_3d_plot

1. 数据生成：与二维散点图类似，通过随机数生成 x、y 和 z 数组的数据。
2. 绘图操作：call pgenv3d 设置三维坐标轴范围，call pgpt3d(ndata, x, y, z, 1) 绘制三维散点图，其中 1 同样表示使用实心圆作为点的符号。

交互式图形绘制

在一些情况下，我们希望能够与绘制的图形进行交互，例如缩放、平移等操作。虽然 PGPLOT 本身的交互功能有限，但可以结合一些其他工具或库来实现更丰富的交互。

结合 X11 事件处理实现简单交互

下面是一个简单示例，展示如何结合 X11 事件处理在 Fortran 绘制的图形上实现交互缩放功能。

program interactive_plot
    use pgplot
    use, intrinsic :: iso_c_binding
    implicit none
    integer :: i, ndata
    real :: x(100), y(100)
    integer(c_int) :: event
    real(c_double) :: x1, y1, x2, y2
    ndata = 100
   ! 生成数据
    do i = 1, ndata
        x(i) = real(i)
        y(i) = sin(x(i) / 10.0)
    end do
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('/XWINDOW')
   ! 设置绘图区域
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv(0.0, real(ndata), -1.2, 1.2, 0, 0)
   ! 绘制线图
    call pgline(ndata, x, y)
   ! 事件循环
    do
        call pgqinq(event, x1, y1, x2, y2)
        if (event == 10) then! 鼠标左键按下事件
            call pgsvp(0.15, 0.95, 0.15, 0.95)
            call pgenv(x1, x2, y1, y2, 0, 0)
            call pgline(ndata, x, y)
        end if
        if (event == 11) exit! 鼠标右键按下事件退出
    end do
   ! 关闭 PGPLOT
    call pgclos()
end program interactive_plot

1. 事件处理：在 do 循环中，通过 call pgqinq(event, x1, y1, x2, y2) 获取 X11 事件。当检测到鼠标左键按下事件（event == 10）时，根据鼠标点击的位置重新设置坐标轴范围，实现缩放功能。当检测到鼠标右键按下事件（event == 11）时，退出事件循环并关闭图形。

图形输出与文件格式

输出到 PostScript 文件

除了在窗口中显示图形，我们还可以将绘制的图形输出到文件，如 PostScript 文件。以下是将二维线图输出到 PostScript 文件的示例。

program plot_to_ps
    use pgplot
    implicit none
    integer :: i, ndata
    real :: x(100), y(100)
    ndata = 100
   ! 生成数据
    do i = 1, ndata
        x(i) = real(i)
        y(i) = sin(x(i) / 10.0)
    end do
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('plot.ps')
   ! 设置绘图区域
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv(0.0, real(ndata), -1.2, 1.2, 0, 0)
   ! 绘制线图
    call pgline(ndata, x, y)
   ! 添加标题和坐标轴标签
    call pglab('X 轴', 'Y 轴', '正弦函数曲线')
   ! 关闭 PGPLOT
    call pgclos()
end program plot_to_ps

在上述代码中，call pgopen('plot.ps') 将输出设备设置为名为 plot.ps 的 PostScript 文件。绘制完成后，会生成该 PostScript 文件，可使用相关的 PostScript 查看器（如 Ghostview）查看。

转换为其他格式

生成的 PostScript 文件可以通过一些工具转换为其他常见的图形格式，如 PNG、JPEG 等。例如，在 Linux 系统下，可以使用 convert 工具（属于 ImageMagick 软件包）将 PostScript 文件转换为 PNG 文件：

convert plot.ps plot.png

这样就可以得到一个更便于在各种平台上查看和分享的 PNG 图像文件。

可视化应用案例

科学数据可视化

在科学研究中，经常需要对实验数据或模拟结果进行可视化。例如，在流体力学模拟中，我们可能得到不同位置的速度数据。假设我们有一个二维网格上的速度数据 u 和 v，可以使用 Fortran 和 PGPLOT 将其可视化。

program fluid_velocity_plot
    use pgplot
    implicit none
    integer :: i, j, nx, ny
    real :: x(100), y(100), u(100, 100), v(100, 100)
    nx = 100
    ny = 100
   ! 生成模拟数据（这里简单示例，实际数据来自模拟）
    do i = 1, nx
        x(i) = (i - 1) / real(nx - 1)
        do j = 1, ny
            y(j) = (j - 1) / real(ny - 1)
            u(i, j) = cos(2 * 3.14159 * x(i)) * sin(2 * 3.14159 * y(j))
            v(i, j) = -sin(2 * 3.14159 * x(i)) * cos(2 * 3.14159 * y(j))
        end do
    end do
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('/XWINDOW')
   ! 设置绘图区域
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv(0.0, 1.0, 0.0, 1.0, 0, 0)
   ! 绘制矢量图
    call pgsls(1)
    call pgsci(1)
    call pgvect(nx, ny, x, y, u, v, 0.02)
   ! 添加标题和坐标轴标签
    call pglab('X 方向', 'Y 方向', '流体速度矢量图')
   ! 关闭 PGPLOT
    call pgclos()
end program fluid_velocity_plot

1. 数据生成：这里简单地通过数学函数生成模拟的速度数据 u 和 v。在实际应用中，这些数据将来自复杂的流体力学模拟程序。
2. 绘图操作：call pgvect(nx, ny, x, y, u, v, 0.02) 用于绘制矢量图，展示流体在不同位置的速度方向和大小。其中 0.02 是矢量箭头的长度缩放因子。

工程数据分析可视化

在工程领域，例如结构力学分析中，我们可能需要对结构的应力分布进行可视化。假设我们有一个二维结构上各点的应力数据 sigma，以下是使用 Fortran 和 PGPLOT 绘制应力分布伪彩色图的示例。

program stress_plot
    use pgplot
    implicit none
    integer :: i, j, nx, ny
    real :: x(200), y(200), sigma(200, 200)
    nx = 200
    ny = 200
   ! 生成模拟数据（实际数据来自结构分析）
    do i = 1, nx
        x(i) = (i - 1) / real(nx - 1)
        do j = 1, ny
            y(j) = (j - 1) / real(ny - 1)
            sigma(i, j) = (x(i) - 0.5) ** 2 + (y(j) - 0.5) ** 2
        end do
    end do
   ! 初始化 PGPLOT
    call pgpap(0.0, 0.0)
    call pgopen('/XWINDOW')
   ! 设置绘图区域
    call pgsvp(0.15, 0.95, 0.15, 0.95)
    call pgenv(0.0, 1.0, 0.0, 1.0, 0, 0)
   ! 绘制伪彩色图
    call pgsci(1)
    call pgcmap(1)
    call pgsurf(nx, ny, x, y, sigma)
   ! 添加颜色条
    call pgscmap(1)
    call pgcbar(0.97, 0.15, 0.99, 0.85)
   ! 添加标题和坐标轴标签
    call pglab('X 方向', 'Y 方向', '应力分布伪彩色图')
   ! 关闭 PGPLOT
    call pgclos()
end program stress_plot

1. 数据生成：这里简单生成了模拟的应力数据。实际应用中，这些数据将来自有限元分析等结构力学计算程序。
2. 绘图操作：call pgsurf 用于绘制表面图，这里根据应力数据 sigma 绘制出伪彩色图。call pgcbar 用于添加颜色条，显示应力值与颜色的对应关系。

通过以上内容，我们详细介绍了在 Fortran 中使用 PGPLOT 库进行图形绘制与可视化的方法，包括二维和三维图形绘制、交互式绘图、图形输出以及在科学和工程领域的应用案例。希望这些内容能帮助你在 Fortran 项目中实现强大的可视化功能。