Hi everyone,
I am trying to plot data from a .nc file with a Lambert native grid, but
in a Mercator projection.
I can't find any working solution for several days: I try first with
rcrm2rgrid function, but when I launch my ncl script, nothing happens
and ncl seems to run for hours..
Here is my script:
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"
load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/shea_util.ncl"
;************************************************
a = addfile(v_file,"r")
t = a->TMP_GDS3_HTGL_10(:,:)
lat2d= a->g3_lat_0
lon2d= a->g3_lon_1
tt = t/10.
dimll = dimsizes(lat2d)
nlat = dimll(0)
nlon = dimll(1)
tt@lat2d = lat2d
tt@lon2d = lon2d
printVarSummary(tt) wks_type = "ps"
wks =gsn_open_wks(wks_type,"tempe")
gsn_define_colormap(wks,"BlWhRe") ; choose color map
res = True ; plot mods desired
res@gsnMaximize = True
res@gsnPaperOrientation = "landscape"
res@gsnBoxMargin = 0.0
res@lbLabelBarOn = False
res@tmXBOn = False ; don't draw axis nor axis labels
res@tmXTOn = False
res@tmYLOn = False
res@tmYROn = False
res@cnLevelSelectionMode = "ManualLevels"
res@cnMinLevelValF = -15.
res@cnMaxLevelValF = 36.
res@cnLevelSpacingF = 1.
res@gsnAddCyclic = False
res@mpDataBaseVersion = "RANGS_GSHHS"
res@mpDataResolution = "FinestResolution"
res@mpPerimOn = False
res@mpProjection = "mercator"
res@mpLimitMode = "LatLon"
res@mpMinLatF = 40
res@mpMaxLatF = 50.3
res@mpMinLonF = -8.02
res@mpMaxLonF = 20.82
lat = new((/nlat/), float)
lon = new((/nlon/), float)
xgrd = rcm2rgrid(lat2d,lon2d,tt,lat,lon,1)
res@gsnAddCyclic = False
plot=gsn_csm_contour_map(wks,xgrd,res)
end
Then I tried to follow the
http://www.ncl.ucar.edu/Applications/Scripts/native_5.ncl script
exemple, but as I want a final mercator projection, so with
one-dimensional arrays that specifies the latitude and longitude
coordinates of the grid instead of two-dimensional in the exemple, I
made some changes:
res@mpProjection = "LambertConformal"
res@mpLimitMode = "Corners"
res@mpLeftCornerLatF = lat2d@corners(0)
res@mpLeftCornerLonF = lon2d@corners(0)
res@mpRightCornerLatF = lat2d@corners(2)
res@mpRightCornerLonF = lon2d@corners(2)
res@mpLambertParallel1F = lat2d@Latin1
res@mpLambertParallel2F = lat2d@Latin2
res@mpLambertMeridianF = lat2d@Lov
res@tfDoNDCOverlay = True
res@gsnDraw = False
plot=gsn_csm_contour_map(wks,tt,res)
getvalues plot
"vpXF" : x
"vpYF" : y
"vpWidthF" : w
"vpHeightF" : h
end getvalues
print("x: " + x + " y: " + y + " w: " + w + " h: " + h)
ndcx = fspan(x, x + w,nlon)
ndcy = fspan(y-h,y,nlat)
eps = 5e-7
ndcx(0) = ndcx(0) + eps
ndcx(nlon-1) = ndcx(nlon -1) - eps
ndcy(0) = ndcy(0) + eps
ndcy(nlat-1) = ndcy(nlat -1) - eps
print( "calculating lat/lon fields")
tndcy = new(nlon, float)
outlat = new((/nlat,nlon/), float)
outlon = new((/nlat,nlon/), float)
do i = 0, nlat-1
tndcy = ndcy(i)
ndctodata(plot,ndcx,tndcy,outlon(i,:),outlat(i,:))
end do
res@mpProjection = "mercator"
res@mpLimitMode = "LatLon"
res@mpMinLatF = 40
res@mpMaxLatF = 50.3
res@mpMinLonF = -8.02
res@mpMaxLonF = 20.82
res@mpMinLonF = i_lon_min
res@mpMaxLonF = i_lon_max
res@tfDoNDCOverlay = False
res@gsnAddCyclic = False
res@sfXArray = outlon(nlat-1,:)
res@sfYArray = outlat(:,nlon-1)
plot=gsn_csm_contour_map(wks,tt,res)
end
But my Postscript output is empty.
The ncdump output of my input file is like this:
netcdf \2010013023 {
dimensions:
g3_x_0 = 339 ;
g3_y_1 = 437 ;
variables:
float g3_lon_1(g3_x_0, g3_y_1) ;
g3_lon_1:La1 = 34.449f ;
g3_lon_1:Lo1 = -10.703f ;
g3_lon_1:Lov = 6.25f ;
g3_lon_1:Dx = 6955.f ;
g3_lon_1:Dy = 7040.f ;
g3_lon_1:Latin1 = 46.8f ;
g3_lon_1:Latin2 = 46.8f ;
g3_lon_1:units = "degrees_east" ;
g3_lon_1:GridType = "Lambert Conformal Secant or
Tangent, Conical or bipolar" ;
g3_lon_1:long_name = "longitude" ;
g3_lon_1:corners = -10.703f, 21.93901f, 29.32028f,
-18.61473f ;
float g3_lat_0(g3_x_0, g3_y_1) ;
g3_lat_0:La1 = 34.449f ;
g3_lat_0:Lo1 = -10.703f ;
g3_lat_0:Lov = 6.25f ;
g3_lat_0:Dx = 6955.f ;
g3_lat_0:Dy = 7040.f ;
g3_lat_0:Latin1 = 46.8f ;
g3_lat_0:Latin2 = 46.8f ;
g3_lat_0:units = "degrees_north" ;
g3_lat_0:GridType = "Lambert Conformal Secant or
Tangent, Conical or bipolar" ;
g3_lat_0:long_name = "latitude" ;
g3_lat_0:corners = 34.449f, 34.66832f, 55.32279f,
55.00012f ;
float g3_rot_2(g3_x_0, g3_y_1) ;
g3_rot_2:note2 = "apply formulas to derive u and v
components relative to earth" ;
g3_rot_2:note1 = "u and v components of vector
quantities are resolved relative to grid" ;
g3_rot_2:formula_v = "Vearth = cos(rot)*Vgrid -
sin(rot)*Ugrid" ;
g3_rot_2:formula_u = "Uearth = sin(rot)*Vgrid +
cos(rot)*Ugrid" ;
g3_rot_2:units = "radians" ;
g3_rot_2:GridType = "Lambert Conformal Secant or
Tangent, Conical or bipolar" ;
g3_rot_2:long_name = "vector rotation angle" ;
float TMP_GDS3_HTGL_10(g3_x_0, g3_y_1) ;
TMP_GDS3_HTGL_10:initial_time = "01/27/2010 (00:00)" ;
TMP_GDS3_HTGL_10:forecast_time_units = "hours" ;
TMP_GDS3_HTGL_10:forecast_time = 28094 ;
TMP_GDS3_HTGL_10:level = 2 ;
TMP_GDS3_HTGL_10:parameter_number = 11 ;
TMP_GDS3_HTGL_10:parameter_table_version = 3 ;
TMP_GDS3_HTGL_10:gds_grid_type = 3 ;
TMP_GDS3_HTGL_10:level_indicator = 105 ;
TMP_GDS3_HTGL_10:coordinates = "g3_lat_0 g3_lon_1" ;
TMP_GDS3_HTGL_10:_FillValue = 1.e+20f ;
TMP_GDS3_HTGL_10:units = "K" ;
TMP_GDS3_HTGL_10:long_name = "Temperature" ;
TMP_GDS3_HTGL_10:center = "French Weather Service -
Toulouse" ;
}
Am I doing something wrong?
Any remark would be very helpful!!
Thanks in advance for your time
Caroline
-- Caroline Serraud - Ingénieur Développement ------------------------------------------------------ METEO CONSULT / La Chaine Météo - Groupe Figaro Département Informatique Domaine de Marsinval F-78540 Vernouillet, FRANCE Tél : 01 39 28 1990 - Fax : 01 39 71 85 31 e-mail : cs@meteoconsult.fr ------------------------------------------------------ Toute la météo sur le Web : http://www.meteoconsult.fr Toute la météo par téléphone : 3201
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Fri Jan 29 02:22:16 2010
This archive was generated by hypermail 2.1.8 : Mon Feb 01 2010 - 08:05:34 MST