计算画图和提取结果出图不一致

lqouc

紫雨微澜 发表于 2013-11-20 18:04

                               
登录/注册后可看大图

难道是fwrite的对精度有影响？但也不会差这么多吧

这个应该不是，图像上看已经不是精度的问题了。

river

紫雨微澜 发表于 2013-11-20 18:05

                               
登录/注册后可看大图

学习了，从来没想过这个问题

那还是对grads的基础知识不熟悉，要用grads计算绘图，就要提供符合grads规定的排列顺序的资料，否则结果就都是错的。可以看一下grads实用手册，上面有说grads要求的数据怎么排列

river

紫雨微澜 发表于 2013-11-20 17:49

                               
登录/注册后可看大图

因为之前是用Z的，但是计算时候pre=lev,用z结果不对，所以就用了判断
'reinit'
'sdfopen I:\daily\air. ...

把经纬度设置换成set x 和set y的形式。还有把你编写的ctl 贴出来

紫雨微澜

river 发表于 2013-11-20 21:36

                               
登录/注册后可看大图

把经纬度设置换成set x 和set y的形式。还有把你编写的ctl 贴出来

dset I:/daily/thse.grd
title mean daily NMC Reanalysis
*dtype grid
*options byteswapped
undef 32766

xdef  144 linear 0 2.5
ydef  73 linear -90 2.5
zdef  8 levels 1000 925 850 700 600 500 400 300
tdef  365 linear 00Z01Jan1979 1440mn
vars 1
thse 8  0 Daily mean thse
endvars

river

紫雨微澜 发表于 2013-11-20 21:44

                               
登录/注册后可看大图

dset I:/daily/thse.grd
title mean daily NMC Reanalysis
*dtype grid

看着也没什么问题，不过缺测值为什么是那个值呢。一般提取资料或者计算输出资料以后缺测值都是grads默认值-9.99e8

紫雨微澜

river 发表于 2013-11-20 22:08

                               
登录/注册后可看大图

看着也没什么问题，不过缺测值为什么是那个值呢。一般提取资料或者计算输出资料以后缺测值都是grads默认值 ...

试过，两个缺测值没有差别

abd

前面人说的那样，garads数据存储读取方式熟悉一下，遇到这种问题的时候建议你自己写个很小的2维或3维数组，每个维度比如用5个数，读写画个简单的图就知道他的规律了。还有，你那个z从1到8那儿是啥意思？有那个必要吗？

紫雨微澜

abd 发表于 2013-11-20 23:07

                               
登录/注册后可看大图

前面人说的那样，garads数据存储读取方式熟悉一下，遇到这种问题的时候建议你自己写个很小的2维或3维数组， ...

当时是想到下面有pres=lev，所以就这么做了，不过又试了下，确实没有必要。但是原来的问题依旧没有解决

'reinit'
'sdfopen I:\daily\air.1979.nc'
'sdfopen I:\daily\rhum.1979.nc'
'set gxout fwrite'
'set fwrite I:\daily\thse3.grd'

    tt=1
  while (tt<=365)
'set t 'tt''

zz=1
  while (zz<=8)
'set z 'zz''
*'set lon 70 140'
*'set lat 15 55'

'define t0=air.1'
'define rh=rhum.2'
'define prs=lev'
'define es=(6.112*exp(17.67*(t0-273.15)/(t0-29.65)))'
'define q=rh*(0.62197*es/(prs-es))/100.'
'define e=prs*q/(0.62197+q)+1e-10'
'define tlcl=55.0+2840.0/(3.5*log(t0)-log(e)-4.805)'
'define theta=t0*pow((1000/prs),(0.2854*(1.0-0.28*q)))'
'define eqt=theta*exp(((3376./tlcl)-2.54)*q*(1.0+0.81*q))'
'define thse=eqt-273.15'
'd thse'

zz=zz+1
endwhile
tt=tt+1
endwhile

'disable fwrite'
'reinit'

yangwawa

{:eb315:}{:eb315:}{:eb315:}{:eb315:}{:eb315:}{:eb315:}{:eb315:}{:eb315:}{:eb315:}

紫雨微澜

yangwawa 发表于 2013-11-21 13:51

                               
登录/注册后可看大图

居然中雨了？表示不可思议啊