在Linux环境下使用OTPS工具箱提取TPXO9海潮模型数据作潮汐水位预报

前期准备

下文中提到的WSL环境，即Windows Subsystem for Linux，具体安装可参考文章 地址

下载OTPS工具箱

这里提供OTPS工具箱的百度网盘下载链接(密码:9z02)：OTPS工具箱
里面包含原始TPXO7.2海潮数据和相关案例。

TPXO9海潮数据请点击下载,原始数据网址 Url。

导入TPXO9数据

我们将下载好的TPXO9数据（共25个）保存在同一个文件夹下，为了区分原始工具包中带有的TPXO7.2海潮数据，将文件夹命名为TPXO9_atlas,我这里是下载了两种数据格式的TPXO9数据（一种是二进制【binary】;一种是nc数据），所以我在TPXO9_atlas的文件夹下再次分了文件夹binary。

我的TPXO9海潮数据保存在 E:\TPXO9_atlas\binary\ 中。

同时有两种方式将TPXO9海潮数据移动到OTPS工具箱中:

一种是直接在Windows界面下直接将TPXO9_atlas文件夹直接移动到OTPS工具箱的DATA文件夹中;

二是通过在WSL中的键入命令命令mv TPXO9_atlas Shell/OTPS/DATA/将其移动到DATA文件夹下。

注意：以上的操作都是要根据自己的实际文件名和文件位置来操作。

提取数据并预报

第一步
在Linux(WSL)中打开到OTPS工具箱下，cd /mnt/e/Shell/OTPS,（以自己实际保存路径为准）

zhuyh@LAPTOP-464DFIKA:OTPS$ ls
COPYRIGHT  Interp  constit.h      extract_HC.f90         extract_local_model.f90  makefile          predict_tide_sample_out  setup.local  weights.h
DATA       README  constit_f90.h  extract_HC_sample_out  lat_lon_time（list）     predict_tide.f90  setup.inp                subs.f90

我们已经在前面导入TPXO海潮模型数据，可以在DATA/文件夹下查看

zhuyh@LAPTOP-464DFIKA:OTPS$ cd DATA/
zhuyh@LAPTOP-464DFIKA:DATA$ ls
Model_Hawaii  Model_tpxo7.2  Model_tpxo9.v1  Model_tpxo9_atlas  TPXO9_atlas  grid_tpxo7.2  h_tpxo7.2  load_file  u_tpxo7.2

第二步
我们可以看到在DATA文件下有一个TPXO9_atlas,里面保存的是TPXO9海潮数据(以自己具体保存路径为准，尽量保存在OTPS/DATA/中，以便于后续的操作)
这时候我们修改控制文件 Model_tpxo9_atlas,键入命令vi Model_tpxo9_atlas:

DATA/TPXO9_atlas/binary/h_*_tpxo9_atlas_30
DATA/TPXO9_atlas/binary/u_*_tpxo9_atlas_30
DATA/TPXO9_atlas/binary/grid_tpxo9_atlas_30
                                                                                                                                                                                                                                                                                                                                                                 

将文件中 h,u,grid 的路径位置修改为自己TPXO9数据的路径位置，保存退出:先 ESC 在键入:wq!

第三步
返回上级目录到 OTPS 下，创建模型运行所需要的经纬度+时间文件，也就是自己需要预报的区域和时间文件，如下所示(逐小时预报，分秒为0即可，以空格为间隔)：

    9.795800   103.471901   2002      3     18      4     29     18
    9.795800   103.471901   2002      3     28      2     28     22
    9.795800   103.471901   2002      6     15     10     15     56
    9.795800   103.471901   2002      6     25      8     15      0
  -29.676001    59.550201   2001      9     18      7     39     51
  -29.676001    59.550201   2001      9     28      5     38     54
    2.021030    72.294998   2002      4     24     11     18     31
    2.021030    72.294998   2002      5      4      9     16     53
    2.021030    72.294998   2002      5     14      7     15     14

下面是我自己写的自动生成 经纬度时间 文件的shell脚本(make_lltime.sh)，可以修改其中的lon,lat范围，时间起始点即可生成想要的数据文件。

#!/bin/bash

# lon
lon_min=117
lon_max=120
# lat
lat_min=10
lat_max=12
# interp
interps=0.1
# time
start_date=20210601
end_date=20210701    # 建议比最后时间多一天

rm ll_time lat_lon
# ==================================================
lons=($(seq ${lon_min} ${interps} ${lon_max}))
lats=($(seq ${lat_min} ${interps} ${lat_max}))

for llat in ${lats[@]}
do
        for llon in ${lons[@]}
        do
                echo "$llat $llon" >> lat_lon
        done
done

while [ $start_date != $end_date ]
do
echo $start_date
for hour in `seq -w 0 23`;do
        awk '{print $1,$2,'${start_date:0:4}','${start_date:4:2}','${start_date:6:2}','${hour}',0,0}' lat_lon >> ll_time
done
let start_date=`date -d "-1 days ago ${start_date}" +%Y%m%d`
done

修改完成后，保存退出！在命令行输入 bash make_lltime,等待脚本执行完成，最后生成的 ll_time 文件就是我们需要的经纬度+时间文件！

第四步

编辑输入文件 setup.inp, vi setup.inp ,在文件的下面有每条指令详细的讲解，这里不多做赘述。

DATA/Model_tpxo7.2         ! 1. tidal model control file
../mt_zhuo/ll_time         ! 2. latitude/longitude/<time> file
z                          ! 3. z/U/V/u/v
                           ! 4. tidal constituents to include
AP                         ! 5. AP/RI
oce                        ! 6. oce/geo
1                          ! 7. 1/0 correct for minor constituents
../mt_zhuo/zeta_mt.out     ! 8. output file (ASCII)

要注意 ：
1.海潮模型控制文件的路径；
2.需要提取数据的经纬度或 需要预报的经纬度+时间 文件
3.需要预报的数据类型，水位(z)，潮流(u)，具体见文件
4.提取或预报所需分潮的名称，所有分潮保持为 空 即可
5,6,7 条可以不用更改，有具体要求可看文件中详解
8.输出文件的路劲位置

根据实际情况修改：

DATA/Model_tpxo9_atlas         ! 1. tidal model control file
./ll_time                ! 2. latitude/longitude/<time> file
z                          ! 3. z/U/V/u/v
                           ! 4. tidal constituents to include
AP                         ! 5. AP/RI
oce                        ! 6. oce/geo
1                          ! 7. 1/0 correct for minor constituents
./zeta_mt.out           ! 8. output file (ASCII)

这里我选取所有分潮进行预报，根据控制文件更改好路径位置后，保存退出！

第五步
现在我们已经有了 经纬度+时间 文件、输入文件 setup.inp
键入命令 make

zhuyh@LAPTOP-464DFIKA:OTPS$ make
gfortran -o predict_tide -fconvert=swap -frecord-marker=4 predict_tide.f90 subs.f90
subs.f90:1404:21:

 1404 |       call mapxy(1,1,dx,dy,dlon,dlat,SLAT,SLON,HEMI)
      |                     1
Warning: Rank mismatch in argument ‘x’ at (1) (rank-2 and scalar) [-Wargument-mismatch]
subs.f90:1416:21:

 1416 |       call mapll(1,1,dlon,dlat,dx,dy,SLAT,SLON,HEMI)
      |                     1
Warning: Rank mismatch in argument ‘lon’ at (1) (rank-2 and scalar) [-Wargument-mismatch]
zhuyh@LAPTOP-464DFIKA:OTPS$ ls
COPYRIGHT  Interp  constit.h      extract_HC.f90         extract_local_model.f90  ll_time   predict_tide      predict_tide_sample_out  setup.local  weights.h
DATA       README  constit_f90.h  extract_HC_sample_out  lat_lon_time（list）     makefile  predict_tide.f90  setup.inp                subs.f90

发现文件中生成了 predict_tide 可执行文件

第六步
键入命令: ./predict_tide < setup.inp，因为TPXO9 的原始数据量较大，需等待片刻！

zhuyh@LAPTOP-464DFIKA:OTPS$ ./predict_tide < setup.inp

 Lat/Lon/Time file:./ll_time
 Predict OCEAN tide
 Interpolate minor constituents

 DATA/TPXO9_atlas/binary/h_m2_tpxo9_atlas_30
 Model:        tpxo9_atlas
 Lat limits:     -90.0166702       90.0166702
 Lon limits:      1.66666675E-02   360.016663
 Constituents: m2  s2  k1  o1  n2  p1  k2  q1  2n2 m4  ms4 mn4
 Predict elevations (m)
 Constituents to include: m2  s2  k1  o1  n2  p1  k2  q1  2n2 m4  ms4 mn4
 Reading model...m2  s2  k1  o1  n2  p1  k2  q1  2n2 m4  ms4 mn4  done
 Results are in ./zeta_mt.out
Note: The following floating-point exceptions are signalling: IEEE_DENORMAL

当出现 Results are in****** 即完成预报计算。
最后键入命令 vi zeta_mt.out 进行查看！