Hi everyone,
Below is the code of multiple typhoon best tracks. However, I can only get all the polymarkers plotted, but not the polylines. I wonder if it is related to the delete function at the end of the code. So how can i further modify the code in order to get the multiple polylines drawn on the same plot? The data used can be seen in the attachment.Thanks.
Regards
Richard
>>
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
begin
diri = "C:/typhoon-best-track-data/2007-wp-full/"
fils = systemfunc ("cd " + diri + " ; ls *.txt")
print(fils)
print(dimsizes(fils))
nfil = dimsizes(fils)
print(dimsizes(nfil))
;fili = "bwp182007.txt"
;********************************
; select sub-regions
;********************************
minlat = 0.0 ;min(lat) ; 0.0 ; deg N ; was 0
maxlat = 50.0 ;max(lat) ; 57.0 ; deg N ; was 57
minlon = 100.0 ; min(lon) -.;-100.0 ; - deg E = deg W ; was -100
maxlon = 165.0 ; max(lon) +10.; -10.0 ; - deg E = deg W ; was -20, then -15
;***** Define some NICE color maps to use.
colors_8 = (/"White","Black","Black","MediumPurple1","MediumPurple3","Blue1",\
"CadetBlue3","Aquamarine2", \
"Gold","Tan1","Sienna1","Tomato","VioletRed1", \
"Yellow","LimeGreen","Grey37","Red","Orange","GoldenRod1", \
"DarkOrange","SteelBlue1","SlateBlue1","SlateGray1", \
"LightSlateBlue","Magenta","DodgerBlue", \
"LightSteelBlue1","Moccasin","LightYellow","LemonChiffon1", \
"CornSilk","LightGoldenrodYellow","Tan","PaleTurquoise3"/)
; Next, set the plotting parameters (which cases to plot)
TypeOn = True ; If you just want the legend to have the basic (H, TS, and TD) designators (* and S), set TypeOn to False
; NOTE -- this will cause the program to grab a different input file -- one in which all the
; extratropical (E), remnant low (L), and wave cases (W) have been taken out; additionally,
; there will not be a distinction between the S and * cases (they will just be plotted according to
; intensity
PlotSubTropicalSegments = True ; If you wish to exclude all subtropical line segments, set to False
; Note: if TypeOn is set to False, the following (W, L, E) cases won't be plotted regardless of their value
PlotWaveSegments = True ; If you wish to exclude all wave segments, set to False
PlotLowSegments = True ; If you wish to exclude all remnant low segments, set to False
PlotExtraTropicalSegments = True ; If you wish to exclude all extratropical line sements, set to False
;********************************
; create plot
;********************************
wks = gsn_open_wks("x11","WP2007")
gsn_define_colormap(wks,colors_8)
res = True
res@gsnMaximize = True
res@gsnDraw = False ; so we can add poly stuff
res@gsnFrame = False ; do not advance frame
; res@gsnPaperOrientation = "portrait"
res@mpDataBaseVersion = "Ncarg4_1" ; Alias 'MediumRes'
res@mpDataSetName = "Earth..1"
; res@mpProjection = "LambertConformal"
res@mpLambertParallel1F = 20.0
res@mpLambertParallel2F = 40.0
res@mpLambertMeridianF = -60.0
res@mpLimitMode = "LatLon"
res@mpMinLatF = minlat
res@mpMaxLatF = maxlat
res@mpMinLonF = minlon
res@mpMaxLonF = maxlon
res@mpFillOn = True ; False to turn off gray continents
res@mpOutlineOn = True ; turn on continental outline
res@mpOutlineBoundarySets = "AllBoundaries"
res@mpLandFillColor = "Tan" ; was "GoldenRod1"
res@mpInlandWaterFillColor = "PaleTurquoise3" ; was "LightBlue1" ; was "PaleTurquoise3"
res@mpOceanFillColor = "PaleTurquoise3" ; was "LightBlue1"
res@mpGeophysicalLineColor = "Grey37"
res@mpGeophysicalLineThicknessF = 0.5
res@mpUSStateLineColor = "Grey37"
res@mpUSStateLineThicknessF = 0.5
res@mpNationalLineColor = "Grey37"
res@mpNationalLineThicknessF = 0.5
res@mpGridAndLimbOn = "True"
res@mpGridAndLimbDrawOrder = "Draw"
res@mpGridMaskMode = "MaskLand"
res@mpGridSpacingF = 5.0
res@mpGridLineColor = "Grey37"
res@tmXBLabelFontHeightF = 0.012 ; 0.005
res@tmXBMajorLengthF = -0.001
res@pmTickMarkDisplayMode = "Always"
;res@tiMainString = ssn + " " + cyclone_name + " (" + year(0) + month(0) + ")"
;***************************************
; plot base map *
;***************************************
plot = gsn_csm_map(wks,res) ; draw one of eight map projections
;********************************
; get data
;********************************
do nf=0,nfil-1,1 ; loop over each file
; do nf=0,2 ; loop over each file
data = asciiread(diri+fils(nf), -1,"string")
printVarSummary(data)
;data = asciiread(diri+fili, -1,"string") ; Info includes type so cases can be weeded out
;tChr = stringtochar(data)
nrows = dimsizes(data)
region = str_get_cols( data, 0, 1)
print(region)
year = str_get_cols( data, 8, 11) ; year,...hour are type "string"
month = str_get_cols( data, 12, 13) ; Returns an array of substrings, given a start and end index
day = str_get_cols( data, 14, 15)
hour = str_get_cols( data, 16, 17)
vmax = stringtoint( str_get_cols( data, 48, 50) )
print (vmax)
lat = stringtofloat(str_get_cols( data, 35, 37))*0.1
;if (str_get_cols( data, 38, 38).eq."S") then
; lat = -lat
;end if
print(lat)
lon = stringtofloat(str_get_cols( data, 41, 44))*0.1
;if (str_get_cols( data, 45, 45).eq."W") then
; lon = -lon
;end if
print(lon)
type = str_get_cols( data, 59, 60) ; "DB", "TD"
name = str_get_cols( data,149,158)
print(type+" "+name)
; Now create arrays that only hold the points for the 00z locations
ind0 = ind(hour.eq."00") ; returns the integer indice if hour = "00"
print(ind0)
print(ind0(0))
if (.not.ismissing(ind0(0))) then ; if ind0(0) don't have a missing value
lat0 = lat(ind0)
lon0 = lon(ind0)
end if
print("lat0="+lat0+" lon0="+lon0)
ind12 = ind(hour.eq."12")
if (.not.ismissing(ind12(0))) then
lat12 = lat(ind12)
lon12 = lon(ind12)
day12 = day(ind12)
end if
print("lat12="+lat12+" lon12="+lon12+" day12="+day12)
;***************************************
; Draw best track history as polylines *
;***************************************
res_poly = True ; polyline mods desired
; create array of dummy graphic variables. This is required, b/c each line
; must be associated with a unique dummy variable.
dum = new(dimsizes(hour),graphic)
ic = 0
do k=0,nrows-2
if(.not.ismissing(lat(k)))
res_poly@gsLineDashPattern = 0
res_poly@gsLineThicknessF = 4
if(.not. TypeOn) then
; Plot hurricane segments
if((vmax(k).gt..64).and.(type(k).eq."*")) then
res_poly@gsLineColor = "Red"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot tropical storm segments
if((vmax(k).lt..64).and.(vmax(k).gt.34).and.(type(k).eq."*")) then
res_poly@gsLineColor = "Yellow"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot tropical depression segments
if((vmax(k).lt..34).and.(type(k).eq."*")) then
res_poly@gsLineColor = "LimeGreen"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Now plot the subtropical segments with dashed lines
; Plot hurricane segments of subtropical systems with dashed red lines
if((vmax(k).gt..64).and.(type(k).eq."S")) then
res_poly@gsLineColor = "Red"
res_poly@gsLineDashPattern = 2
res_poly@gsLineDashSegLenF = 0.1
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot tropical storm segments of subtropical systems with dashed yellow lines
if((vmax(k).lt..64).and.(vmax(k).gt.34).and.(type(k).eq."S")) then
res_poly@gsLineColor = "Yellow"
res_poly@gsLineDashPattern = 2
res_poly@gsLineDashSegLenF = 0.1
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot tropical depression segments of subtropical systems with dashed dark LimeGreen lines
if((vmax(k).lt..34).and.(type(k).eq."S")) then
res_poly@gsLineColor = "LimeGreen"
res_poly@gsLineDashPattern = 2
res_poly@gsLineDashSegLenF = 0.1
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
else ; if(TypeOn) -- plot pure tropical segments as solid lines (red for hurricane, etc.).
; Plot super typhoon segments
if((vmax(k).ge..64).and.(type(k).eq."ST")) then
res_poly@gsLineColor = "Purple"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot hurricane segments
if((vmax(k).ge..64).and.(type(k).eq."TY")) then
res_poly@gsLineColor = "Red"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot tropical storm segments
if((vmax(k).lt..64).and.(type(k).eq."TS").and.(vmax(k).ge.34)) then
res_poly@gsLineColor = "Orange"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot tropical depression segments
if((vmax(k).lt..34).and.(type(k).eq."TD")) then
res_poly@gsLineColor = "Yellow"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot tropical disturbance segments
if((type(k).eq.."DB").and.(vmax(k).lt.34)) then
res_poly@gsLineColor = "Green"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
if(PlotWaveSegments) then
; Plot tropical wave segments
if(type(k).eq."WV") then
res_poly@gsLineColor = "LimeGreen"
res_poly@gsLineDashPattern = 2
res_poly@gsLineDashSegLenF = 0.1
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
end if
if(PlotLowSegments) then
; Plot remanant low segments
if(type(k).eq."LO") then
res_poly@gsLineColor = "SlateBlue"
res_poly@gsLineDashPattern = 2
res_poly@gsLineDashSegLenF = 0.1
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
end if
if(PlotSubTropicalSegments) then
; Plot subtropical storm segments
if((vmax(k).gt..34).and.(type(k).eq."SS")) then
res_poly@gsLineColor = "DarkOrange"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
; Plot subtropical depression segments
if((vmax(k).lt..34).and.(type(k).eq."SD")) then
res_poly@gsLineColor = "Blue1"
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
end if
if(PlotExtraTropicalSegments) then
; Plot extratropical segments
if(type(k).eq."EX") then
res_poly@gsLineColor = "Black"
res_poly@gsLineDashPattern = 2
res_poly@gsLineDashSegLenF = 0.1
dum(ic) = gsn_add_polyline(wks,plot,lon(k:k+1),lat(k:k+1),res_poly)
ic = ic + 1
end if
end if
end if
end if
end do
;***************************************
; Draw black dots at lat/lon itersections
;***************************************
;***************************************
; Draw black polymarkers at 00Z locations and white polymarkers at 12z locations
;***************************************
res_mark = True
res_mark@gsMarkerIndex = 1 ; polymarker style
res_mark@gsMarkerSizeF = 0.02 ; polymarker size - was 0.012
res_mark@gsMarkerColor = "Black" ; change marker color
duma = new(1,graphic)
duma = gsn_add_polymarker(wks,plot,lon0,lat0,res_mark)
; Now draw white polymarkers at 12z locations
res_mark@gsMarkerIndex = 1
res_mark@gsMarkerSizeF = 0.02 ; polymarker size - was 0.012
res_mark@gsMarkerColor = "White" ; change marker color
dumb = new(1,graphic)
dumb = gsn_add_polymarker(wks,plot,lon12,lat12,res_mark)
;***************************
; Plot some text labels *
;***************************
txres = True
; Now draw the day labels at each 00Z marker
txres@txFontHeightF = 0.015 ; 0.008
txres@txJust = "TopCenter"
dumz = new(nrows,graphic)
n12 = dimsizes(day12)
;;do k=0,nrows-1
do k=0,n12-1
dumz(k) = gsn_add_text(wks,plot,day12(k),lon12(k),lat12(k)-0.25,txres)
end do
; Draw storm number (actually these are strings) at beginning of track
txres@txFontHeightF = 0.02 ; 0.0044
txres@txJust = "CenterLeft"
txres@txPerimOn = True
txres@txPerimColor = "Black"
txres@txPerimThicknessF = 1.0
txres@txPerimSpaceF = 0.4
txres@txBackgroundFillColor = "White"
; gsn_text(wks,plot,ssn,lon(0)+0.6,lat(0),txres)
; Now draw them at the end
txres@txJust = "BottomRight"
; gsn_text(wks,plot,ssn,lon(nrows-1),lat(nrows-1),txres)
if (.not. TypeOn) then
txres@txFontHeightF = 0.02 ; 0.006
txres@txJust = "CenterCenter"
txres@txPerimSpaceF = 0.7
txres@txPerimThicknessF = 1.5
txres@txConstantSpacingF = 1.0
txres@txFontAspectF = 1.5
txres@txFontThicknessF = 0.8
gsn_text_ndc(wks,"~F02~ Remnant low, tropical wave, and ~C~ extratropical stages are not shown. ~C~~C~ Subtropical stages are dashed.",0.715,0.680,txres)
end if
txres@txFontHeightF = 0.005
txres@txPerimThicknessF = 1.5
txres@txPerimSpaceF = 1.0
txres@txConstantSpacingF = 1.0
txres@txFontAspectF = 1.5
txres@txFontThicknessF = 0.8
;******************
; Plot a legend *
;******************
lgres = False
lgres@lgAutoManage = False
lgres@vpWidthF = 0.2 ; was 0.1 ; width of legend (NDC)
lgres@vpHeightF = 0..2 ; was 0.1 ; height of legend (NDC)
; lgres@lgBottomMarginF = 0.17 ; was 0.25
lgres@lgPerimFill = 0 ; Use solid fill (0) instead of the default hollow fill
lgres@lgPerimFillColor = "Background"
; lgres@lgBoxMajorExtentF = 0.4
lgres@lgBoxMinorExtentF = 0.2 ; controls how wide the box holding the legend items (lines and markers) can be in relation to legend
; lgres@lgBoxBackground = "PaleTurquoise3"
lgres@lgMonoItemType = False ; indicates that we wish to set the item types individually
lgres@lgMonoMarkerIndex = False
lgres@lgMonoLineThickness = False
lgres@lgMonoMarkerThickness = False
lgres@lgMonoMarkerSize = False
; position fine elements of legend relative to some positional values:
xlegend = 0.66 ; 0.820
ylegend = 0.35 ; 0.290
if(TypeOn) then
lgres@lgLabelFont = 0
lgres@lgLabelFontHeightF = 0.012 ; 0.004
lgres@lgLabelFontAspectF = 1.2
lgres@lgLabelConstantSpacingF = 0.0
lgres@lgItemCount = 11
lgres@lgLineDashSegLenF = 0.1
lgres@lgItemTypes = (/"Markers","Markers","Markers","Lines","Lines","Lines","Lines","Lines","Lines","Lines","Lines"/)
lgres@lgMarkerIndexes = (/ 1, 4, 16, 0, 0, 0, 0, 0, 0, 0, 0/)
lgres@lgLineThicknesses = (/ 0.1, 0.1, 0.1, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0, 4.0/)
lgres@lgMarkerColors = (/ "White", "Black", "Black"/)
lgres@lgMarkerSizes = (/ 0.0001, 0.004, 0.004/)
lgres@lgLineColors = (/ "White", "Black", "Black",\
"SlateBlue1", "LimeGreen", "Black", "Blue1", "Dark Orange", "LimeGreen", "Yellow", "Red" /) ; colors for legend lines
lgres@lgDashIndexes = (/ 0, 0, 0,\
2, 2, 2, 0, 0, 0, 0, 0 /) ; dash indexes
; gsn_legend_ndc(wks,10, (/"Tropical Cyclone No.","1200 UTC Position/Date","0000 UTC Position",\
; "Remnant Low","Tropical Wave","Extratropical","Subtropical Dep.","Subtropical Storm","Tropical Dep.","Tropical Storm (T)","Hurricane (H)"/),xlegend,ylegend,lgres)
yoffset_day = -0.0858
yoffset_storm = -0.0910
else
lgres@lgLabelFont = 0
lgres@lgLabelFontHeightF = 0.004
lgres@lgLabelFontAspectF = 1.2
lgres@lgLabelConstantSpacingF = 0.0
lgres@lgItemCount = 6
lgres@lgItemTypes = (/"Markers","Markers","Markers","Lines" ,"Lines" ,"Lines"/)
lgres@lgLineColors = (/"White", "Black" ,"Black" ,"LimeGreen","Yellow","Red" /) ; colors for legend lines
lgres@lgLineThicknesses = (/ 0.1, 0.1, 0.1, 4.0, 4.0, 4.0/)
lgres@lgDashIndexes = (/ 0, 0, 0, 0, 0, 0/)
lgres@lgMarkerIndexes = (/ 0, 4, 16, 0, 0, 0/)
lgres@lgMarkerColors = (/"White", "Black" , "Black"/)
lgres@lgMarkerSizes = (/ 0.0001, 0.004, 0.004/)
; gsn_legend_ndc(wks,5,(/"Tropical Cyclone No.","1200 UTC Position/Date","0000 UTC Position",\
; "Tropical Dep.","Tropical Storm (T)","Hurricane (H)"/),xlegend,ylegend,lgres)
yoffset_day = -0.0750 ; offsets for the case where we don't plot all the extra cases
yoffset_storm = -0.0875
end if
; Now draw a day label on the legend
txres@txPerimOn = False
txres@txFontHeightF = 0.01
txres@txPerimSpaceF = 0.0
txres@txJust = "CenterCenter"
; gsn_text_ndc(wks,"21",xlegend + 0.0183, ylegend + yoffset_day,txres)
; Draw storm number on the legend
txres@txFontHeightF = 0.004
txres@txPerimOn = True
txres@txPerimColor = "Black"
txres@txPerimThicknessF = 1.0
txres@txPerimSpaceF = 0.4
txres@txBackgroundFillColor = "White"
; gsn_text_ndc(wks,"3",xlegend + 0.0142, ylegend + yoffset_storm,txres)
; drawNDCGrid(wks) ; Use to put on a grid for fine-tuning placement
delete(region)
delete(data) ; delete all arrays that may change size
delete(year)
delete(month)
delete(day)
delete(hour)
delete(vmax)
delete(lat)
delete(lon)
delete(type)
delete(name)
delete(ind0)
delete(lat0)
delete(lon0)
delete(lat12)
delete(lon12)
delete(day12)
delete(ind12)
delete(dum)
delete(dumz)
end do
draw(wks)
frame(wks)
end
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
This archive was generated by hypermail 2.1.8 : Thu Dec 03 2009 - 09:53:22 MST