用ncl画Geogrid之后的landuse 图

arthurmaple

跑WRS的时候，做了一点改动，采用了modis+3种城市细分的静态数据。
也就是说绝大多数土地利用都应该是modis20类，其中拓展了三种新的类型，在WPS里面定义为31，32，33

想对geogrid出来的图形进行ncl可视化

我希望呈现的是，31类为brown颜色，32类为red颜色，33类为black 颜色
但是不知道为什么最后32类也被显示成黑色了。。。而且颜色表和字体重合很严重

请问各位同学，关于这个奇怪的现象有什么建议吗？小弟初学ncl，暂时没有头绪，不知道要做怎么处理才好（可以确定的是，中心区域LU_INDEX一定是32，我用panoply看过了）
怀疑是因为lu读取的数据不是连续的（1~20之后就是31,32,33）带来的问题，但是我也不知道。。。。还望大神不吝指点

代码如下：
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/csm/contributed.ncl"   

begin

  f0=addfile("/home/mode/WPS/geo_em.d02.nc","r")
  lon=f0->XLONG_M(0,:,:)
  lat=f0->XLAT_M(0,:,:)
  lu=f0->LU_INDEX(0,:,:)
  lu@description = ""
  lu@units = ""
  printVarSummary(lon)
  printVarSummary(lat)
  printVarSummary(lu)

  info     = (/ " 1  Evergreen Needleleaf", \  ; n=0
                " 2  Evergreen Broadleaf ", \
                " 3  Deciduous Needleleaf", \
                " 4  Deciduous Broadleaf ", \
                " 5  Mixed Forest        ", \
                " 6  Closed Shrublands   ", \
                " 7  Open Shrublands     ", \
                " 8  Woody Savannas      ", \
                " 9  Savannas            ", \
                "10  Grasslands          ", \
                "11  Permanent Wetlands  ", \
                "12  Croplands           ", \
                "13  Urban and Built-up  ", \
                "14  Cropland Mosaics    ", \
                "15  Snow and Ice        ", \
                "16  Bare Soil and Rocks ", \
                "17  Water Bodies        ", \
                "31  low indensity residential area", \
                "32  high indensity residential area", \
                "33  commercial area"  /) ; n=19

ninfo = dimsizes(info)

colors     = (/ "white","black","darkgreen", \  ; 1
"darkolivegreen", \      ;2
"darkkhaki", \     ;3
"green2", \ ;4
"limegreen", \    ;5
"darkolivegreen1", \         ;6
"greenyellow", \                            ;7
"darkolivegreen3", \    ;8
"darkolivegreen3", \                ;9
"darkolivegreen2", \         ; 10
"orange1", \            ;11
"coral2", \        ;12
"red", \  ;13
"chocolate1", \ ;  14
"lightyellow1", \  ;15
"gold1", \ ; 16
"blue", \ ; 17
"white", \   ;31
"red", \   ;32
"black" /) ; 33

;************************************************
; create plot
;************************************************

  wks = gsn_open_wks("ps","anduse")  


gsn_define_colormap(wks,colors)     

  res                  = True               ; plot mods desired
  res@gsnDraw          = False
  res@gsnFrame         = False
  res@gsnMaximize      = True
  res@gsnAddCyclic     = False              ; Don't add a cyclic point, important-zhiyong
  res@cnFillOn         = True               ; color Fill
  res@cnFillMode       = "RasterFill"       ; Raster Mode
  res@cnLinesOn        = False              ; Turn off contour lines

  res@cnLevelSelectionMode = "ExplicitLevels" ; set explict contour levels
  res@cnLevels         = integertobyte( ispan(2,ninfo,1) )

  res@lbLabelPosition  = "Center"           ; label position
  res@lbLabelAlignment = "BoxCenters"       ; label orientation
  res@lbLabelStrings   = ispan(1,ninfo,1)  

  res@pmLabelBarOrthogonalPosF = -0.08    ; move up smidge
  res@pmTickMarkDisplayMode = "Always"

  res@tmXBMode = "Automatic"
  res@tmYLMode = "Automatic"
  res@tmXTLabelsOn = False
  res@tmYRLabelsOn = False
  res@tmXTMajorLineColor  = 0
  res@tmYRMajorLineColor  = 0


  res@mpCenterLonF     = 0                  ; set map center
  res@mpFillOn         = False
  res@mpOutlineOn      = False              ; set default boundaries
  res@tfDoNDCOverlay   = True
  res@mpProjection      = "LambertConformal"
  res@mpLambertParallel1F = 30.0
  res@mpLambertParallel2F = 60.0
  res@mpLambertMeridianF  = 113.0        ; central meridian
  res@mpLimitMode         = "LatLon"

  res@sfXArray =lon
  res@sfYArray =lat

  res@mpMinLatF =  min(lat)
  res@mpMaxLatF =  max(lat)
  res@mpMinLonF =  min(lon)
  res@mpMaxLonF =  max(lon)

  res@tmXBMode = "Automatic"
  res@tmYLMode = "Automatic"
  res@tmXBLabelFontHeightF  = 0.01

;  plot = gsn_csm_contour_map_ce(wks, lu, res) ; create plot
   plot = gsn_csm_contour_map(wks, lu, res)

;
;add text

  rtxt = True
  rtxt@txJust        = "CenterLeft"
  rtxt@txFontHeightF = 0.010
                            ; 3 rows x 6 columns of text
  n  = -1
  xx = 0.05                 ; arbitrary
  do ncol=0,5
     yy = 0.22
    do nrow=0,2
       n = n+1
       gsn_text_ndc (wks,info(n),xx,yy,rtxt)
       yy = yy - 3*rtxt@txFontHeightF
    end do
     xx = xx + 0.155
  end do

  draw(plot)
  frame(wks)

arthurmaple

有大神能帮忙看看吗

arthurmaple

解决了，确实是怀疑的那个问题
通过将info里面33种全部填满（不能17之后就直接31 32 33）
以及颜色表那里也全部填满（18~30我都设置了白色）

目前可以正确出图

唯一的问题是，底下的色表和现实的LU对应文字会多出不需要的18~30这么多

目前还不知道怎么解决，我打算图用这个图，但是文字说明自己P上去吧。。。

南方是有个我呀

"跑WPS的时候，做了一点改动，采用了modis+3种城市细分的静态数据。
也就是说绝大多数土地利用都应该是modis20类，其中拓展了三种新的类型，在WPS里面定义为31，32，33"

请问您如何在WPS里定义新类型31\32\33的呢?

arthurmaple

南方是有个我呀 发表于 2019-5-28 20:32
"跑WPS的时候，做了一点改动，采用了modis+3种城市细分的静态数据。
也就是说绝大多数土地利用都应该是mod ...

这个说起来不那么容易，你搜一下替换下垫面，有很多类似的文章和帖子

jiajun102481

想请教下增加3种城市细分数据，看有些帖子说除了修改namelist.input外，还要强制设置process_domian_module.F中的num_land_cat？不知道思路对不对？

shenyi

arthurmaple 发表于 2019-6-21 18:11
这个说起来不那么容易，你搜一下替换下垫面，有很多类似的文章和帖子

你好，请问你修改土地分类数据，都修改了哪些部分啊，是修改了geog.吗

shenyi

请问下，修改土地利用时，原13类修改为31-33，那原13类你保留了吗

arthurmaple

shenyi 发表于 2020-1-4 18:12
请问下，修改土地利用时，原13类修改为31-33，那原13类你保留了吗

我保留了一部分，因为精细下垫面资料不全，有的只好保留

shenyi

arthurmaple 发表于 2020-1-6 10:39
我保留了一部分，因为精细下垫面资料不全，有的只好保留

能否这样理解，要研究的城市部分，将13进行修改为31-33，该城市之外的部分，13可以不修改