[Thread Prev][Thread Next][Index]

Re: [ferret_users] contour/fill with logarithmic scaled Y - Axis



Forgot to cc the list:

Hi Anja,

You need to re-grid you data onto an axis that consists of the logarithmic values of your z-coordinates. Assuming that you just want the log values of the sigma coordinates (rather than re-gridding to a 4D true pressure field), something like this should work:

yes? define axis/depth/from_data/name=lnpz/z ln(z[gz=MY_DATA])
yes? let x_lnp=MY_DATA[gz=lnpz@asn]

where MY_DATA could be your temperature trend. I have this in a short script (p2lnp.jnl, see below), where $1 replaces MY_DATA. To re-grid to true pressure, you need to use the ZAXREPLACE function - search the email archives for some examples of this.

The plotting requires a little more tinkering. The two scripts below (one for plotting, one for labelling the new y axes) are based on something from Glenn Carver (U. Cambridge). The latter script is good for most ranges that you may want to specify (e.g. 1000-10 hPa), but you can always pull out the "label" and "aline" commands and run them individually to give particular axis definitions (e.g. if you want to label 763 hPa for some reason).

Hope this helps,

Paul

SCRIPTS:

(a) Plotting
---------------
! plot_lnp.jnl
!
! Description: Plot a z axis plot with ln(p)
!
! DETAILS:
! Plot z axis plot with log(p) as vertical coordinate. Need to run one of 
! conversion scripts (theta2lnp or p2lnp) first to give an appropriately
! gridded variable. Use in conjunction with lnp_labels to get the z axis
! labelled correctly. 
!
! ARGUMENTS:
!    $1    Command to use to plot (e.g. fill/lev=(-10,10,1) )
!    $2    Variable to plot
!    $3    zmin (e.g. 1000 hPa) - for plot range
!    $4    zmax (e.g.   10 hPa) - for plot range 
!    $5    Title of the plot
!
! EXAMPLE:
! Plot zonal mean temperature from model output on sigma hybrid pressure
! coordinates.
!
! yes? use temperature.nc
!       yes? go p2lnp temperature                 !Get log P z axis; creates x_lnp variable
!       yes? go plot_lnp fill/l=10/lev=(270,330,10) x_lnp "Temperature"
!       yes? ppl xlab "Latitude"
!       yes? ppl fill
!       yes? go lnp_labels 1000 10 0.11          !Do the labelling
!

! Define input variables
let u_lnp = $2
let/quiet zmin=$3
let/quiet zmax=$4

! Set the position or turn off the labels
! axlabp labx, laby  -1 = bottom/left, 0 = no label, 1 = top/right
ppl axlabp -1, 0

! Suppress large and small tics on y axis. Leave x as default
ppl tics 0.125, 0.25, 0, 0

! Do the plot
set variable/title="$5" u_lnp
set region/z=`ln(zmin)`:`ln(zmax)`
$1/set u_lnp

(b) Y-axis Labelling
---------------------------
! lnp_labels.jnl
!
! Description: Label z axis for log P plot
!
! DETAILS:
! Control labelling for z axis on log P type plot. Use in conjunction
! with plot_lnp.jnl. 
!
! ARGUMENTS:
!    $1      zmin (e.g. 1000 hPa)
!    $2      zmax (e.g.   10 hPa)
!    $3      size for z axis labels (e.g. 0.11)
!

let/quiet zmin = $1
let/quiet zmax = $2
let/quiet size = $3

! Do the major tic marks on the z axis (log)

! First set the grid for the plot so we can use proper 'users' coords
! in the repeat command.
let/quiet logmin = `log(zmin)`
let/quiet logmax = `log(zmax)`

define axis/x=`logmax`:`logmin`:1 xx
define grid/x=xx gg
set grid gg

! Do the tic for the powers of 10 first
repeat/x=`logmax`:`logmin`:1 ppl aline 1,($ppl$xmin),`ln(10^x)`,-92,`ln(10^x)`

! Do the tics for the values '2', '3', '5' and '7' inbetween the major powers of 10
! We don't do all the tics as it'll get too crowded.
let zxmin = logmin-1

! The zmin value may not be an exact power of ten in which case we might need to
! do a 2*10^x value. Should really test to see if this exceeds the $ppl$ymin
! before plotting
!
!(a) Left hand axis
if `int(logmin) ne logmin AND ln(2*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 ppl aline 1,($ppl$xmin),`ln(2*10^x)`,-92,`ln(2*10^x)`
else 
repeat/x=`logmax`:`int(logmin)-1`:1 ppl aline 1,($ppl$xmin),`ln(2*10^x)`,-92,`ln(2*10^x)`
endif
if `int(logmin) ne logmin AND ln(3*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 ppl aline 1,($ppl$xmin),`ln(3*10^x)`,-92,`ln(3*10^x)`
else
repeat/x=`logmax`:`int(logmin)-1`:1 ppl aline 1,($ppl$xmin),`ln(3*10^x)`,-92,`ln(3*10^x)`
endif
if `int(logmin) ne logmin AND ln(3*10^int(logmin)) le ($ppl$ymax)` then
repeat/x=`logmax`:`int(logmin)`:1 ppl aline 1,($ppl$xmin),`ln(5*10^x)`,-92,`ln(5*10^x)`
else
repeat/x=`logmax`:`int(logmin)-1`:1 ppl aline 1,($ppl$xmin),`ln(5*10^x)`,-92,`ln(5*10^x)`
endif
if `int(logmin) ne logmin AND ln(5*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 ppl aline 1,($ppl$xmin),`ln(7*10^x)`,-92,`ln(7*10^x)`
else
repeat/x=`logmax`:`int(logmin)-1`:1 ppl aline 1,($ppl$xmin),`ln(7*10^x)`,-92,`ln(7*10^x)`
endif

!(b) Right hand axis - just tics
if `int(logmin) ne logmin AND ln(2*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 ppl aline 1,($ppl$xmax),`ln(2*10^x)`,92,`ln(2*10^x)`
else 
repeat/x=`logmax`:`int(logmin)-1`:1 ppl aline 1,($ppl$xmax),`ln(2*10^x)`,92,`ln(2*10^x)`
endif
if `int(logmin) ne logmin AND ln(3*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 ppl aline 1,($ppl$xmax),`ln(3*10^x)`,92,`ln(3*10^x)`
else
repeat/x=`logmax`:`int(logmin)-1`:1 ppl aline 1,($ppl$xmax),`ln(3*10^x)`,92,`ln(3*10^x)`
endif
if `int(logmin) ne logmin AND ln(3*10^int(logmin)) le ($ppl$ymax)` then
repeat/x=`logmax`:`int(logmin)`:1 ppl aline 1,($ppl$xmax),`ln(5*10^x)`,92,`ln(5*10^x)`
else
repeat/x=`logmax`:`int(logmin)-1`:1 ppl aline 1,($ppl$xmax),`ln(5*10^x)`,92,`ln(5*10^x)`
endif
if `int(logmin) ne logmin AND ln(5*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 ppl aline 1,($ppl$xmax),`ln(7*10^x)`,92,`ln(7*10^x)`
else
repeat/x=`logmax`:`int(logmin)-1`:1 ppl aline 1,($ppl$xmax),`ln(7*10^x)`,92,`ln(7*10^x)`
endif

! Do the labelling
let/quiet rightjust = 1

! powers of 10 first
repeat/x=`logmax`:`logmin`:1 label -93, `ln(10^x)`, `rightjust`, 0, `size`, "`10^x`"

if `int(logmin) ne logmin AND ln(2*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 label -93,`ln(2*10^x)`,`rightjust`,0,`size`, "`2*10^x`"
else 
repeat/x=`logmax`:`int(logmin)-1`:1 label -93,`ln(2*10^x)`,`rightjust`,0,`size`, "`2*10^x`"
endif
if `int(logmin) ne logmin AND ln(3*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 label -93,`ln(3*10^x)`,`rightjust`,0,`size`, "`3*10^x`"
else
repeat/x=`logmax`:`int(logmin)-1`:1 label -93,`ln(3*10^x)`,`rightjust`,0,`size`, "`3*10^x`"
endif
if `int(logmin) ne logmin AND ln(3*10^int(logmin)) le ($ppl$ymax)` then
repeat/x=`logmax`:`int(logmin)`:1 label -93,`ln(5*10^x)`,`rightjust`,0,`size`, "`5*10^x`"
else
repeat/x=`logmax`:`int(logmin)-1`:1 label -93,`ln(5*10^x)`,`rightjust`,0,`size`, "`5*10^x`"
endif
if `int(logmin) ne logmin AND ln(5*10^int(logmin)) le ($ppl$ymax)` then 
repeat/x=`logmax`:`int(logmin)`:1 label -93,`ln(7*10^x)`,`rightjust`,0,`size`, "`7*10^x`"
else
repeat/x=`logmax`:`int(logmin)-1`:1 label -93,`ln(7*10^x)`,`rightjust`,0,`size`, "`7*10^x`"
endif


! Label y axis (not sure why ferret doesn't do this?)
let/quiet posy = (ln(zmin) + ln(zmax)) / 2
let/quiet centre = 0
let/quiet angle = 90
!label -108,`posy`,`centre`, `angle`, `size*2`, "@ASPressure (hPa)"     !Sometimes too big- PY
label -108,`posy`,`centre`, `angle`, `size*1.3`, "@ASPressure (hPa)"    !Better

  
On Apr 15, 2009, at 4:54 AM, ajna.West@xxxxxx wrote:

Dear ferreters,

I use the following script to produce a contour plot. (see attachment)

How can I convert the Y - Axis into a logarithmic scaled Y - Axis?

Thanks in advance for any help,

Anja  


My script:

USE 1960_2004_ALLzmb.nc

fill/nolab var130
contour/over/nolab var130

LABEL -105,50000,0,90,.18 Z (level)
LABEL 0,109500,0,0,.18 Latitude
LABEL 0,-10500,0,0,.18 Temperaturtrend (K)


This is the grid of the netCDF file:

yes? show grid
Default grid for DEFINE VARIABLE is ABSTRACT
Last successful data access was on grid GQI1
   GRID GQI1
name       axis              # pts   start                end
LON       LONGITUDE            1mr   0E                   0E
LAT       LATITUDE            48 i   87.159S              87.159N
LEV       Z (level)           22 i-  1000                 100000
TIME      T (day as %Y%m       1 r   541102               541102
<Temperaturtrend.gif>


-----
Paul Young

Chemistry and Climate Processes
Chemical Sciences Division
NOAA Earth System Research Laboratory
325 Broadway R/CSD8
Boulder CO 80305-3328
USA

Tel:   +1 303-497-4711
Fax:   +1 303-497-5686



-----
Paul Young

Chemistry and Climate Processes
Chemical Sciences Division
NOAA Earth System Research Laboratory
325 Broadway R/CSD8
Boulder CO 80305-3328
USA

Tel:   +1 303-497-4711
Fax:   +1 303-497-5686


[Thread Prev][Thread Next][Index]

Contact Us
Dept of Commerce / NOAA / OAR / PMEL / TMAP

Privacy Policy | Disclaimer | Accessibility Statement