AWX格式tbb在micaps，NCL和Meteoinfo的数值问题

zeldazhe

研究了一下用NCL读取和画AWX格式的tbb数据。现在图像以及叠加在正确经纬度信息上基本解决。但是发现数值大小问题无法解决。
读取出来的tbb数据值在-127-128之间。
参考：运用NCL处理风云卫星2E的AWX格式数据的总结
http://bbs.06climate.com/forum.php?mod=viewthread&tid=48371&fromuid=29971
把tbb数据做了如下处理：
do i=0,511
          do j=0,511
            index=tbb(i,j)
            if (index .lt. 0) then
                tbb(i,j) = tbb(i,j)+255
            else
                tbb(i,j) = tbb(i,j)
            end if
          end do
并参考其色标进行画图：
res@cnLevelSelectionMode                    = "ExplicitLevels"
        ; *************************************色标设置****************************
            colors = (/(/255,255,255/),(/  0,  0,  0/),(/33,33,33/),(/68,65,68/),(/132,128,134/),(/201,192,201/),(/230,222,232/),\
                    (/235,166,195/),(/239,142,160/),(/244,118,123/),(/245,73,73/),(/220,79,69/),(/153,118,78/),(/181,145,80/),\
                    (/224,187,89/),(/239,202,92/),(/219,207,130/),(/205,212,165/),(/191,217,199/),(/186,218,209/),(/163,213,234/),\
                    (/156,203,228/),(/146,190,220/),(/144,187,217/),(/140,182,214/),(/132,172,208/),(/130,169,206/),(/116,149,193/),\
                    (/108,139,186/),(/103,132,182/),(/103,132,182/),(/95,121,175/),(/95,121,175/), (/84,106,165/),(/74,93,157/),\
                    (/68,86,152/),(/65,82,149/),(/59,73,143/),(/49,60,133/),(/45,55,129/), (/43,52,127/),(/36,45,121/),(/34,43,120/),\
                    (/35,42,118/),(/33,39,116/),(/27,32,110/),(/22,25,104/),(/19,21,101/),(/17,18,98/),(/13,13,94/),(/11,11,92/),\
                    (/10,11,91/),(/10,11,91/),(/10,11,91/),(/10,11,91/),(/10,10,91/),(/11,12,92/),(/12,12,93/),(/12,12,93/),\
                    (/11,11,92/),(/10,10,91/),(/10,10,91/)/)*1.0 ; we multiply by 1 to make colors float      

            res@cnLevels                                = fspan(50,250,59)
            res@cnFillColors                            = ispan(61,2,1)
            res@cnConstFLabelTextDirection              = "Down"
            cmap = colors/255.             ; normalize (required by NCL)
图片和在micaps显示色调的基本一致。但是数值不一致
处理后数据值在0-255之间。同时，此数据的高值区对应的是micaps读图同样颜色的低值区。

meteoinfo也是一样的结果，数值在0-255之间。这个值该如何转换为摄氏度或者开尔文温度？
以下三幅图分别是micaps，ncl和meteoinfo的图，注意色标对应的值。

micaps

NCL

meteoinfo

zeldazhe

附上我使用的awx格式文件 ANI_IR1_R01_20140626_1530_FY2D.AWX (1.38 MB, 下载次数: 5)

2017-6-16 16:14 上传

点击文件名下载附件

59093558@qq.com

请教一下， 您的脚本中 1200*1200格点， 对应所用到的land.dat 经纬度对照表 是如何定义的？ 万分感谢帮忙 🙏

59093558@qq.com

弱弱的问下， 请问能不能学习下您的ncl脚本呢？

59093558@qq.com

load "./shapefile_utils.ncl"
begin
system("date")

if (.not. isvar("awxFilePath")) then
    awxFilePath = "/Users/nannan/Desktop/SATELLITE/FY2E/L1/IR1/EQUAL/201712/ANI_IR1_R04_20171226_2330_FY2E.AWX"
end if
if (.not. isvar("lenRecord")) then
    lenRecord = 1900   ;记录长度
end if
if (.not. isvar("numData")) then
    numData = 1300   ;产品数据占用记录数
end if
if (.not. isvar("numHead")) then
    numHead = 3   ;文件头占用记录数
end if
if (.not. isvar("lonCellDistance")) then
    lonCellDistance = 0.04997368421052631   ;经度间隔
end if
if (.not. isvar("latCellDistance")) then
    latCellDistance = 0.04994615384615384   ;纬度间隔
end if
if (.not. isvar("lonStart")) then
    lonStart = 50.06997368421053   ;起始经度
end if
if (.not. isvar("latStart")) then
    latStart = -4.910053846153846   ;起始纬度
end if
if (.not. isvar("imgWidth")) then
    imgWidth = 1900   ;图像宽度
end if
if (.not. isvar("imgHeight")) then
    imgHeight = 1300   ;图像高度
end if

if (.not. isvar("outputFile")) then      ; is outputFile command line?
  outputFile = "/Users/nannan/Desktop/99";
end if

f=fbindirread(awxFilePath,0,(/numData+numHead,lenRecord/),"byte")
dataArray = byte2flt(f(numHead:(numData+numHead-1),:))

data=new((/numData,lenRecord/),float)
do i=0,numData-1
    index = numData-i-1
    do j=0,lenRecord-1
        temp = dataArray(i,j)
        if (temp .lt. 0) then
            temp = temp+255
        end if
        data(index,j) = temp
    end do
end do   

;***********************************设置数据属性****************************
olon = new(lenRecord,float)
olat = new(numData,float)

do i=0, numData-1
    olat(i)= latStart + latCellDistance * i
end do

do j=0, lenRecord-1
  olon(j)=lonStart + lonCellDistance * j
end do

rectMinLat = olat(0)
rectMaxLat = olat(numData-1)
rectMinLon = olon(0)
rectMaxLon = olon(lenRecord-1)

olon!0 = "lon";
olon@long_name = "lon";
olon@units = "degrees_east";
olon&lon = olon;

olat!0 = "lat"
olat@long_name = "lat";
olat@units = "degrees_north";
olat&lat=olat;

data!0        = "lat"
data!1         = "lon"
data&lon = olon
data&lat = olat
data@_FillValue            = -9999.9

;*************************************设置工作台****************************
res                                         = True
res@gsnDraw                                 = False
res@gsnFrame                                = False
res@gsnMaximize                             = True
;res@gsnPaperOrientation                     = "portrait"
res@gsnAddCyclic                            = False                                 
res@pmTickMarkDisplayMode                   = "Always"
;*************************************网格设置****************************
res@mpGridAndLimbOn                     = False
res@mpLimitMode                         = "LatLon"
res@mpMinLatF                           =   rectMinLat        
res@mpMaxLatF                           =  rectMaxLat
res@mpMinLonF                           =  rectMinLon
res@mpMaxLonF                           = rectMaxLon

res@trGridType                              = "TriangularMesh"
res@cnFillOn                                = True
res@cnFillMode                              = "RasterFill"
res@cnLinesOn                               = False
res@cnLineLabelsOn                          = False
res@cnLevelSelectionMode                    = "ExplicitLevels"

res@lbLabelBarOn                     = False
; *************************************色标设置****************************
colors = (/(/255,255,255/),(/  0,  0,  0/),(/33,33,33/),(/68,65,68/),(/132,128,134/),(/201,192,201/),(/230,222,232/),\
          (/235,166,195/),(/239,142,160/),(/244,118,123/),(/245,73,73/),(/220,79,69/),(/153,118,78/),(/181,145,80/),\
          (/224,187,89/),(/239,202,92/),(/219,207,130/),(/205,212,165/),(/191,217,199/),(/186,218,209/),(/163,213,234/),\
          (/156,203,228/),(/146,190,220/),(/144,187,217/),(/140,182,214/),(/132,172,208/),(/130,169,206/),(/116,149,193/),\
          (/108,139,186/),(/103,132,182/),(/103,132,182/),(/95,121,175/),(/95,121,175/), (/84,106,165/),(/74,93,157/),\
          (/68,86,152/),(/65,82,149/),(/59,73,143/),(/49,60,133/),(/45,55,129/), (/43,52,127/),(/36,45,121/),(/34,43,120/),\
          (/35,42,118/),(/33,39,116/),(/27,32,110/),(/22,25,104/),(/19,21,101/),(/17,18,98/),(/13,13,94/),(/11,11,92/),\
          (/10,11,91/),(/10,11,91/),(/10,11,91/),(/10,11,91/),(/10,10,91/),(/11,12,92/),(/12,12,93/),(/12,12,93/),\
          (/11,11,92/),(/10,10,91/),(/10,10,91/)/)*1.0 ; we multiply by 1 to make colors float      

res@cnLevels                                = fspan(50,250,59)
res@cnFillColors                            = ispan(61,2,1)

cmap = colors/255.             ; normalize (required by NCL)
;*************************************绘制设置****************************
;创建工作台
wks_type = "png"
wks_type@wkWidth                                = imgWidth          ;工作台宽度
wks_type@wkHeight                               = imgHeight        ;工作台高度
wks = gsn_open_wks(wks_type,outputFile)

contour = gsn_csm_contour_map_ce(wks,data,res)
gsn_define_colormap(wks,cmap)

ynShp="shp/prov.shp"
lnres        = True  
lnres@minlat = 20.99  
lnres@maxlat = 29.54  
lnres@minlon = 97.35  
lnres@maxlon = 106.26  
lnres@gsLineColor      = "blue"
lnres@gsLineThicknessF = 2.0   
shpid1 = gsn_add_shapefile_polylines(wks,contour,ynShp,lnres)

draw(contour)
frame(wks)
delete([/wks,res,contour/])
system("date")
end

Sozo

你好，这个问题解决了吗？

zeldazhe

Sozo 发表于 2018-3-30 15:23
你好，这个问题解决了吗？

没有，后来放弃了，下圆盘HDF格式的数据了

zeldazhe

Sozo 发表于 2018-3-30 15:23
你好，这个问题解决了吗？

没有，后来放弃了，下圆盘HDF格式的数据了