Re: gradient computation

Hi Dennis,

I think making the relationship 1deg -> 111.123e3 m for all latitudes could
result in errors away from 45 deg lat.
I made some tests comparing the meridional gradient of temperature with two
sets of data, both from reanalysis. Please, see below:

; read in data
f = addfile( "air.1979.nc", "r" ) ; data with dlat
fixed
;f = addfile( "air.2m.gauss.1979.nc", "r" ) ; data with dlat not fixed
T = short2flt( f->air(0,0,:,:) ) ; temperature at
lowest level and for the first time
lat = f->lat
lon = f->lon

; deg to rad
g2r = acos(-1.)/180

; for data with dlon and dlat fixed, equal to 2.5 degrees
; for data with dlat not fixed, the latitude difference used can introduces
errors
; first way to compute the gradient => directly by center_finite_diff
dT_dy = center_finite_diff( T(lon|:,lat|:), (lat(2)-lat(1))*g2r, False, 0 )
* 6.37e-6

; second way to compute gradient => computing variations of temperature and
latitude in y direction
dTy = center_finite_diff( T(lon|:,lat|:), 1., False, 0 ) ; variation
of temperature in y direction
dy = center_finite_diff( lat, 1., False, 0 ) *g2r ;
variation of latitude in y direction (in radians)
dyND = conform( T(lon|:,lat|:), dy, 1 ) ; same
dimensions of dTy
dTy_dyND = (dTy / dyND) * 6.37e-6 ; the
GRADIENT

; copy coodinates
copy_VarCoords( T(lon|:,lat|:), dTy_dyND )
copy_VarCoords( T(lon|:,lat|:), dT_dy )

Ploting 'dT_dy' and 'dTy_dyND' should result in the same graphics, since
dlat = 2.5, it's fixed.

But, when dlat isn't fixed, the graphics will have very little differences
(this differences can be seen by ploting individual grid point values), but
I think if the gradient will be used in other computations, the error can be
propagated. So, for non-fixed dlat, the second way to compute gradient will
be better.

Best regards,

Mateus

2009/7/6 Dennis Shea <shea_at_ucar.edu>

> I must leave but this might be more clear
>
> ; compute K/deg_of_lat
> dT_dLat = center_finite_diff( T( lon|:, lat|: ), lat, False, 0 )
>
> ; convert to K/m: one deg lat @ 45N is approx 111.123e3 meters
> dT_dLat = dT_dLat/111.123e3
> dT_dLat_at_units = "K/m"
>
> Of course, this could be done in one statement.
>
> If this works, please post to the ncl-talk question.
>
> THX
> D
> >
>
> Mateus Teixeira wrote:
>
>> Dennis,
>>
>> Couldn't it be as below
>>
>> dy = center_finite_diff( lat, 1., False, 0 ) * 6.37e6 * acos( -1. )/180.
>> dyND = conform( T, dy, 1 )
>> dTy = center_finite_diff( T( lon|:, lat|: ), 1., False, 0 )
>>
>> dTy_dyND = dTy / dyND
>>
>>
>> Mateus
>>
>>
>>
>> 2009/7/6 Dennis Shea <shea_at_ucar.edu <mailto:shea_at_ucar.edu>>
>>
>> What I sent you was incorrect.
>> I will have to think about this.
>>
>> D
>>
>> Mateus Teixeira wrote:
>>
>> Dear NCL users,
>>
>> I have a temperature field that is not equally spaced in y, but
>> it is in x. I'm not sure how to compute its gradient in y
>> direction. I'm wondering if I can apply this
>>
>> dphi = center_finite_diff( lat, 1., False, 0 ) * 6.37e6 * acos(
>> -1. )/180. dT_dy = center_finite_diff( T( lon|:, lat|: ),
>> dphi, False, 0 )
>>
>> or I just use
>>
>> dT_dy = center_finite_diff( T( lon|:, lat|: ), lat *
>> acos(-1.)/180., False, 0 ) / 6.37e6
>>
>> Thanks.
>>
>>
>> Best regards,
>>
>> -- Mateus da Silva Teixeira
>>
>> Registered Linux User #466740
>>
>>
>>
>> ------------------------------------------------------------------------
>>
>> _______________________________________________
>> ncl-talk mailing list
>> List instructions, subscriber options, unsubscribe:
>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>
>>
>>
>>
>> --
>> Mateus da Silva Teixeira
>>
>> Registered Linux User #466740
>>
>

-- 
Mateus da Silva Teixeira
Registered Linux User #466740

_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Mon Jul 06 2009 - 20:37:49 MDT