Matlab TRMM_3B42数据格式转换（nc转tif，附代码）

[2021年4月21日更新]
TRMM_3B42_Daily数据【TRMM (TMPA) Precipitation L3 1 day 0.25 degree x 0.25 degree V7 (TRMM_3B42_Daily)】无需转置；
TRMM_3B42RT_Daily数据【TRMM (TMPA-RT) Near Real-Time Precipitation L3 1 day 0.25 degree x 0.25 degree V7 (TRMM_3B42RT_Daily)】无需转置；
TRMM_3B42RT_3h数据【TRMM (TMPA-RT) Near Real-Time Precipitation L3 3 hour 0.25 degree x 0.25 degree V7 (TRMM_3B42RT)】需要转置。
TRMM_3B42_3h的数据【TRMM (TMPA) Rainfall Estimate L3 3 hour 0.25 degree x 0.25 degree V7 (TRMM_3B42)】为HDF格式，本文只考虑nc转tif格式，暂不考虑HDF格式。
三种数据的比较见文末。

[2021年7月14日更新]
更新代码中的Refference = georasterref('RasterSize', size(preciSets), 'Lonlim', [-180 180], 'Latlim', [-50 50])语句，代码会根据数据的经纬度范围自适应设置。
另外，TMPA-RT的纬度范围为[-60 60]。

原代码：

Refference = georasterref('RasterSize', size(preciSets), 'Lonlim', [-180 180], 'Latlim', [-50 50]);

更新代码：

Refference = georasterref('RasterSize', size(preciSets), ...
        'Lonlim', double([floor(min(lon(:))) ceil(max(lon(:)))]), ...
        'Latlim', double([floor(min(lat(:))) ceil(max(lat(:)))]) ...
    );

---

TRMM 3B42数据转为tif格式。TRMM数据的下载步骤请参考博客—TRMM降水数据下载步骤。

---

在nc文件转tif文件时，常会遇到影像旋转的问题，遇到此类问题常需要使用rot90()函数将数据进行转置。但是并不是所有的nc数据都需要进行转置，在此之前需要判断数据的行列数以及经纬度（lon和lat）的长度。

【TRMM_3B42RT_3h数据的转置不能使用rot90()函数，直接用B = A'即可】

在读取nc文件的数据（data）部分以及经纬度后进行判断。如果data部分的数组为A x B，而lat的长度为A，lon的长度为B时，此时不需要进行转置，例如TRMM 3B42的nc4格式的数据中，data部分为400 x 1440，lat的长度为400，lon的长度为1440，此时无需转置；反之则需要转置。

另一个问题是通过georasterref()函数设置地理栅格数据参考对象（类）时，需要注意'Lonlim'和'Latlim'的取值范围，此处的设置需要参考lon和lat。例如lon的范围为[-179.875 179.875]，lat的范围为[-49.875 49.875]时，‘Lonlim’ 和 ‘Latlim’ 的设置应该为'Lonlim', [-180 180], 'Latlim', [-50 50]。

TRMM_3B42RT_3h数据转为tif格式数据的参考代码以及结果如下：

%%
clc;
clear;

%%
ncpath = 'D:\20200115\3B42RT.2013092200.7.nc4';
ncinf = ncinfo(ncpath);
preciSets = ncread(ncpath, 'precipitation');
preciSets = preciSets';
lon = ncread(ncpath, 'lon');
lat = ncread(ncpath, 'lat');
OutputPath = 'D:\20200115\3B42RT.2013092200.7.tif';
%地理栅格数据参考对象（类）
Refference = georasterref('RasterSize', size(preciSets), ...
        'Lonlim', double([floor(min(lon(:))) ceil(max(lon(:)))]), ...
        'Latlim', double([floor(min(lat(:))) ceil(max(lat(:)))]) ...
    );
%写成GeoTif格式%写成GeoTif格式
geotiffwrite(OutputPath, preciSets, Refference);   

---

TRMM_3B42_Daily数据 / TRMM_3B42RT_Daily数据转为tif格式数据的参考代码以及结果如下：

%%
clc;
clear;

%%
ncpath = 'D:\20200115\3B42_Daily.20000823.7.nc4';
ncinf = ncinfo(ncpath);
preciSets = ncread(ncpath, 'precipitation');
lon = ncread(ncpath, 'lon');
lat = ncread(ncpath, 'lat');
OutputPath = 'D:\20200115\3B42_Daily.20000823.7.tif';
%地理栅格数据参考对象（类）
Refference = georasterref('RasterSize', size(preciSets), ...
        'Lonlim', double([floor(min(lon(:))) ceil(max(lon(:)))]), ...
        'Latlim', double([floor(min(lat(:))) ceil(max(lat(:)))]) ...
    );
%写成GeoTif格式%写成GeoTif格式
geotiffwrite(OutputPath, preciSets, Refference);   

---

补充对比：

展示的是2013年9月22日的数据，参考图为https://journals.ametsoc.org/view/journals/wcas/11/2/wcas-d-18-0053_1.xml的Figure 1(a)

---

TRMM_3B42_Daily数据，单位为mm d^-1，DN值范围为[0, 332.7]：

---

TRMM_3B42RT_Daily数据，单位为mm d^-1，DN值范围为[0, 228.96]：

---

TRMM_3B42RT_3h（单位为mm hr^-1）合成的Daily数据（3B42RT.2013092200.7.nc4 / 3B42RT.2013092203.7.nc4 / 3B42RT.2013092206.7.nc4 / 3B42RT.2013092209.7.nc4 / 3B42RT.2013092212.7.nc4 / 3B42RT.2013092215.7.nc4 / 3B42RT.2013092218.7.nc4 / 3B42RT.2013092221.7.nc4），DN值范围为[-19.5, 92.49]（乘以3之后即为daily）：

另外， TRMM_3B42RT_3h数据存在负值。
3B42RT.2013092203.7.nc4表示2013年9月22日3时±1.5小时内的降水，详情请参考TRMM_3B42RT_3h数据主页。