Hi,
This script is hard to read, as in my mailer, the indentation
reflecting code structure is all messed up. However, I *think* what
is going on is that all the polyline drawing code surrounded by an
if-clause of the form:
"if (.not.TypeOn) then"...
but prior to that, TypeOn is set to True, so the polyline drawing code
would never execute.
Rick
On Tue, 1 Dec 2009 11:14:23 +0800 (HKT)
Li Richard <lcy114@yahoo.com.hk> wrote:
> 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
>
>
>
> Yahoo!香港提供網上安全攻略，教你如何防範黑客! 請前往
>http://hk.promo.yahoo.com/security/ 了解更多!
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Mon Nov 30 20:51:38 2009
This archive was generated by hypermail 2.1.8 : Thu Dec 03 2009 - 09:53:22 MST