ncl在底图上画闪电

黄小仙儿

想在设定的经纬度区域内画闪电图，闪电资料有三列，纬度，经度，强度，强度大于0就在相应经纬度位置上画+，小于零就画-，总是实现不了，把脚本贴出来求大家指导。


load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_code.ncl"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/gsn_csm.ncl"
;load "$NCARG_ROOT/lib/ncarg/nclscripts/wrf/WRFUserARW.ncl"
begin
  type = "pdf"
  wks = gsn_open_wks(type,"lightning16")
; gsn_define_colormap(wks,"precip_11lev")

; a = addfile("/home/Huanglei/data/d032"+".nc","r")      
  res = True
  res@MainTitle = "REAL-TIME WRF"

  mpres  = True  ; Map resources
  mpres@mpOutlineOn = False  ; Turn off map outlines
  mpres@mpFillOn    = False  ; Turn off map fill
  mpres@mpGridAndLimbOn = True
  pltres = True ; Plot resources
  pltres@PanelPlot  = True   ; Tells wrf_map_overlays not to remove overlays
  mpres@pmTickMarkDisplayMode = "Always"
  lon1=117.80
  lat1=24.20
  lon2=118.35
  lat2=24.4
;---Zoom in on area of interest
mpres@mpLimitMode           = "Corners"
  mpres@mpLeftCornerLatF    = lat1
  mpres@mpLeftCornerLonF    = lon1
mpres@mpRightCornerLatF   = lat2
mpres@mpRightCornerLonF   = lon2
;---Create map.
  map = gsn_csm_map(wks,res2)
  return(map)


;>============================================================<
;                      add China map
;>------------------------------------------------------------<
     
  shp_name1    = "/home/Huanglei/map/China/diquJie_polyline.shp"
  lnres                  = True
  lnres@gsLineColor      = "gray25"
  lnres@gsLineThicknessF = 0.5   
id = gsn_add_shapefile_polylines(wks,map,shp_name1,lnres)
  shp_name2    = "/home/Huanglei/map/China/cnmap/cnhimap.shp"
  prres=True
  prres@gsLineThicknessF = 2.0      
  prres@gsLineColor = "black"
  plotcn3 = gsn_add_shapefile_polylines(wks,map,shp_name2,prres)
; draw(plot)       ; This will draw the map and the shapefile outlines.
      
  ascii_filename = "/home/Huanglei/data/new16jia.txt"
  seismic = asciiread(ascii_filename,(/521,3/),"float")   
  
    y = seismic(:,0)  ; Column 1 of file contains X values.
    x = seismic(:,1)  ; Column 2 of file contains Y values.
    z = seismic(:,2)  ; Column 3 of file contains Z values.
  txres2  = True
  txres2@txFont  = 0.01
  txres2@txFontHeightF =0.01
  txres2@txFontColor = "Red"
idx = ind(z .gt. 0)
print(idx)
if .not. all(ismissing(idx))
    str = new(dimsizes(idx), "string")
    str = "+"
    txdum1 = gsn_add_text(wks, map, str, x(idx),y(idx), txres2)
end if
txres2@txFontColor = "Blue"
idx := ind(z .lt. 0)
if .not. all(ismissing(idx))
    str := new(dimsizes(idx), "string")
    str = "-"   
    txdum2 = gsn_add_text(wks, map, str, x(idx),y(idx), txres2)
end if



  draw(plot)
   frame(wks)
   ;delete(opts)
;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;;
         
end

guoke199006

感谢楼主分享

黄小仙儿

guoke199006 发表于 2015-8-7 11:07
感谢楼主分享

我这个是错的，放上来让大家指错的