Vectors

[Patrick Brockmann <brock@xxxxxxxxxxxxxxxxxx>]

Hi Ferreters,

I haven't seen answer to a mail posted about plotting vectors
and problems related to the use of the vector command with
4 arguments (vector u,v,long_grid,lat_grid).

To resume

yes? use monthly_navy_winds.cdf
yes? vector UWND[l=1],VWND[l=1]

works properly but

yes? vector
UWND[l=1],VWND[l=1],x[g=UWND]*y[g=UWND]*0+x[g=UWND],x[g=UWND]*y[g=UWND]*0+y[g=UWND]

not! So impossible for me to produce curvilinear (and evenmore
scattered) vectors plots.

To solve this, you could use the attached script which can be adapted
to produce other projections.

Hope this will help
Regards
Patrick

--
Patrick BROCKMANN (LSCE/CNRS-CEA Saclay)
Laboratoire des Sciences du Climat et de l'Environnement
Unite mixte de recherche CNRS-CEA
L'Orme des Merisiers - Bat.701 - CEA/Saclay
91191 GIF sur YVETTE cedex - FRANCE
Tel : 01 69 08 86 79 - Fax : 01 69 08 77 16
mailto:brock@lsce.saclay.cea.fr

\cancel mode verify

!**************************************************************
! Description: plot over vectors 
!
! Usage: go plot_vectors u v lon lat [skip] [fill]
!                  
!
! Notes:  default skip = 1 (all vectors)
!         u = zonal wind component
!         v = meridional wind component
!         lon = starting longitude of the vector 
!         lat = starting latitude of the vector 
!         Those arguments can be either 2D or 1D
!        
!         u, v, lon, lat can be either vectors or arrays
!
!         A vector key could be added but the scale is the 
!         same as the previous plot.
!
! Example: use monthly_navy_winds.cdf
!          go basemap " " " " " " " "  hollow
!          go plot_vectors UWND[l=1] VWND[l=1] x[g=UWND]*y[g=UWND]*0+x[g=UWND] x[g=UWND]*y[g=UWND]*0+y[g=UWND]
!
!          go basemap " " " " " " " "  hollow
!          go plot_vectors UWND[l=1] VWND[l=1] x[g=UWND]*y[g=UWND]*0+x[g=UWND] x[g=UWND]*y[g=UWND]*0+y[g=UWND] 10
!
!          go basemap " " " " " " " "  hollow
!          go plot_vectors UWND[l=1] VWND[l=1] x[g=UWND]*y[g=UWND]*0+x[g=UWND] x[g=UWND]*y[g=UWND]*0+y[g=UWND] 10 fill
!
!          go basemap " " " " " " " "  hollow
!          let XPOS = {50,130,240,300,140,90,200,210,110,70}
!          let YPOS = {10,-50,40,20,-10,0,40,70,-35,50}
!          let U = {10,3,-12,8,15,6,-8,9,17,-8}
!          let V = {-3,5,7,2,-3,10,-15,11,9,-12}
!          go plot_vectors U V XPOS YPOS " " fill 
!
! Calls:
!
! Author: Patrick Brockmann
! Contact: brock@lsce.saclay.cea.fr
!
! Creation: October 2000
! Version: 1.0
! History:
! Modification:
!
!**************************************************************

query/ignore $1%<Use: go plot_vectors u v lon lat [skip] [fill]%
query/ignore $2%<Use: go plot_vectors u v lon lat [skip] [fill]% 
query/ignore $3%<Use: go plot_vectors u v lon lat [skip] [fill]% 
query/ignore $4%<Use: go plot_vectors u v lon lat [skip] [fill]% 

if `$5%1% lt 1` then
        message/continue "Error\! Must have option skip > 0"
        exit
endif

!**************************************************************
def sym PV_PI = 3.14159265
def sym PV_DEG2RAD = `($PV_PI)/180`

!**************************************************************
let PV_u=xsequence($1)
let PV_v=xsequence($2)
let PV_lon=xsequence($3)
let PV_lat=xsequence($4)

def sym PV_nbvectors=`PV_u,return=iend`
!say `($PV_nbvectors)`

!-----------------------
! Select 1 vector on $5 (default = all)
let PV_select0=mod(i[i=1:($PV_nbvectors)]-1,`$5%1%`) + 1 
let PV_select1=if PV_select0 eq 1 then 1

let PV_u1 = PV_u * PV_select1
let PV_v1 = PV_v * PV_select1
let PV_lon1 = PV_lon * PV_select1
let PV_lat1 = PV_lat * PV_select1

!-----------------------
! take into account ratio (= ratioworld * ratioviewp)
! to respect in everye case the head angle choosen 
!-----------------------
def sym PV_ratioworld = `(($YAXIS_MAX)-($YAXIS_MIN))/(($XAXIS_MAX)-($XAXIS_MIN))`
def sym PV_ratioviewp = `($PPL$XLEN)/($PPL$YLEN)`
def sym PV_ratio = `($PV_ratioworld)*($PV_ratioviewp)`
!say ($PV_ratioworld) ($PV_ratioviewp) ($PV_ratio) 

!-----------------------
! vector points
!           e 
!         / | \
!        a -|- b
!           | 
!           | 
!           | 
!           s 
!-----------------------
def sym PV_head_size = 0.3                        ! in percentage
def sym PV_head_angle = `30*($PV_DEG2RAD)`            ! in degrees
let PV_size = ($PV_head_size) * (PV_u1^2 + (PV_v1/($PV_ratio))^2)^.5
let PV_alpha = atan2(PV_v1/($PV_ratio),PV_u1)

!-- point a
let PV_headlon_a0 = cos(PV_alpha - ($PV_head_angle)) * PV_size * 1/cos(($PV_head_angle))
let PV_headlat_a0 = sin(PV_alpha - ($PV_head_angle)) * PV_size * 1/cos(($PV_head_angle))

let PV_headlon_a1 = PV_lon1 + PV_u1 - PV_headlon_a0  
let PV_headlat_a1 = PV_lat1 + PV_v1 - PV_headlat_a0 * ($PV_ratio) 

!-- point b
let PV_headlon_b0 = sin(PV_alpha + ($PV_head_angle) - ($PV_PI)/2) * PV_size * 1/cos(($PV_head_angle))
let PV_headlat_b0 = cos(PV_alpha + ($PV_head_angle) - ($PV_PI)/2) * PV_size * 1/cos(($PV_head_angle))

let PV_headlon_b1 = PV_lon1 + PV_u1 + PV_headlon_b0 
let PV_headlat_b1 = PV_lat1 + PV_v1 - PV_headlat_b0 * ($PV_ratio)

!-- make the polygon vectors 
let PV_dum = PV_u*0 + j[j=1:6]

!-----------------------
! normal vectors or filled head vectors
if $6%false|fill>true|*>false% then

  ! head filled in black
  !-- gather points (pt_s--pt_e--pt_a--pt_b--pt_e--pt_s) 
  let PV_headlon0 = if PV_dum eq 1 then PV_lon1 else PV_dum
  let PV_headlon1 = if PV_dum eq 2 then PV_lon1+PV_u1 else PV_headlon0
  let PV_headlon2 = if PV_dum eq 3 then PV_headlon_a1 else PV_headlon1
  let PV_headlon3 = if PV_dum eq 4 then PV_headlon_b1 else PV_headlon2
  let PV_headlon4 = if PV_dum eq 5 then PV_lon1+PV_u1 else PV_headlon3
  let PV_headlon5 = if PV_dum eq 6 then PV_lon1 else PV_headlon4

  let PV_headlat0 = if PV_dum eq 1 then PV_lat1 else PV_dum
  let PV_headlat1 = if PV_dum eq 2 then PV_lat1+PV_v1 else PV_headlat0
  let PV_headlat2 = if PV_dum eq 3 then PV_headlat_a1 else PV_headlat1
  let PV_headlat3 = if PV_dum eq 4 then PV_headlat_b1 else PV_headlat2
  let PV_headlat4 = if PV_dum eq 5 then PV_lat1+PV_v1 else PV_headlat3
  let PV_headlat5 = if PV_dum eq 6 then PV_lat1 else PV_headlat4

  polygon/nolab/over/pal=black/line/fill PV_headlon5, PV_headlat5, i[i=1:($PV_nbvectors)] 

else

  ! normal vectors
  !-- gather points (pt_s--pt_e--pt_a--pt_e--pt_b--pt_e) 
  let PV_headlon0 = if PV_dum eq 1 then PV_lon1 else PV_dum
  let PV_headlon1 = if PV_dum eq 2 then PV_lon1+PV_u1 else PV_headlon0
  let PV_headlon2 = if PV_dum eq 3 then PV_headlon_a1 else PV_headlon1
  let PV_headlon3 = if PV_dum eq 4 then PV_lon1+PV_u1 else PV_headlon2
  let PV_headlon4 = if PV_dum eq 5 then PV_headlon_b1 else PV_headlon3
  let PV_headlon5 = if PV_dum eq 6 then PV_lon1+PV_u1 else PV_headlon4

  let PV_headlat0 = if PV_dum eq 1 then PV_lat1 else PV_dum
  let PV_headlat1 = if PV_dum eq 2 then PV_lat1+PV_v1 else PV_headlat0
  let PV_headlat2 = if PV_dum eq 3 then PV_headlat_a1 else PV_headlat1
  let PV_headlat3 = if PV_dum eq 4 then PV_lat1+PV_v1 else PV_headlat2
  let PV_headlat4 = if PV_dum eq 5 then PV_headlat_b1 else PV_headlat3
  let PV_headlat5 = if PV_dum eq 6 then PV_lat1+PV_v1 else PV_headlat4

  polygon/nolab/over/line PV_headlon5, PV_headlat5, i[i=1:($PV_nbvectors)] 

endif
!-----------------------

!**************************************************************
cancel sym PV_*
cancel var PV_* 

!**************************************************************
set mode/last verify