Re: Data from a Lambert native grid reprojected on a Mercator grid.

From: Dennis Shea <shea_at_nyahnyahspammersnyahnyah>
Date: Fri Jan 29 2010 - 11:56:25 MST

Re:

    lat = new((/nlat/), float)
    lon = new((/nlon/), float)
    xgrd = rcm2rgrid(lat2d,lon2d,tt,lat,lon,1)

You have defined lat & lon to be length nlat and nlon.
However, there are no lat or lon values!!!

Th function

        function rcm2rgrid (
                lat2d [*][*] : numeric,
                lon2d [*][*] : numeric,
                fi : numeric,
                lat [*] : numeric,
                lon [*] : numeric,
                Option : numeric
        )

say that lat an lon are one dimensional arrays
that conatain the lat & lon the data should be interpolated to.

eg:
         lat = fspan(20, 70, nlat)
         lat@units = "degrees_north"

         lon = fspan(40, 230, nlon)
         lon@units = "degrees_east"

Bon chance!

D

Caroline Serraud wrote:
> 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
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk
Received on Fri Jan 29 11:56:30 2010

This archive was generated by hypermail 2.1.8 : Mon Feb 01 2010 - 08:05:34 MST