[ncl-talk] Question on interpolation
Andrew Kren - NOAA Affiliate
andrew.kren at noaa.gov
Wed Sep 4 10:04:13 MDT 2019
Thanks Dennis! That worked perfectly! You may close this ticket.
Thanks,
On Tue, Sep 3, 2019 at 6:21 PM Dennis Shea <shea at ucar.edu> wrote:
> Sorry if what I wrote was not clear.
>
> A sketch of a script follows.
> ===========================
> dirHWRF = "..."
> filHWRF = "..."
> pthHWRF = dirHWRF + filHWRF
> nfil = *dimsizes*(filHWRF)
>
> ; Trick to have NCL do book-keeping of 'time'.
> ; Have NCL's addfiles create a 'time' array that spans all times
> ; It is assumed that there is one time step per file
> ;;setfileoption("grb","SingleElementDimensions",(/"Initial_time",
> "Forecast_time"/))
> *setfileoption*
> <https://www.ncl.ucar.edu/Document/Functions/Built-in/setfileoption.shtml>("grb","SingleElementDimensions",
> "Initial_time")
>
> fadd = *addfiles*
> <https://www.ncl.ucar.edu/Document/Functions/Built-in/addfiles.shtml>
> (pthHWRF,"r")
> time = fadd*[:]*->initial_time0_hours ; times acriss all files
> time!0 = "time"
> ntim = *dimsizes*(time)
> ymdh = c*d_calendar*(time,-3)
> ymdh!0 = "time"
> print(time+" "+ymdh)
> print("------------")
>
> if (ntim.ne.nfil) then
> print("nim="+ntim+" : nfil="+nfil: should be equal)
> exit
> end if
>
> plev = *ff[0]*->lv_ISBL0 ; get from 1st file [invariant]
> plev!0 = "plev"
> print(plev)
> delete(fadd)
>
> ;---Sonde locations
>
> SLAT = (/lat1, lat2, lat3/) ; one or more Sonde lat
> SLON = (/lon1, lon2, lon3/)
> npts = dimsizes(SLAT)
>
> xsonde = new((/ntim,klev,npts/),"float", -1e10) ; allocate space for
> variable
> ;;ysonde = new((/ntim,klev,npts/),"float", -1e10)
> ;;zsonde = new((/ntim,klev,npts/),"float", -1e10)
>
> do nf=0,nfil-1 one time step per file; lat_0, lon_0 coordinates may
> change
> f = *addfile*(pthHWRF*(nf)*,"r")
> x = f->X ; X( initial_time0_hours0, lv_ISBL0, lat_0, lon_0 )
> xsonde(nf,:,:)) = *linint2_points_Wrap*
> <https://www.ncl.ucar.edu/Document/Functions/Contributed/linint2_points_Wrap.shtml>(x&lon_0,
> x&lat_0, x, False,SLAT,SLON 0)
> end do
>
> xsonde!0 = "time"
> xsonde&time = time
> xsonde!1 = "plev"
> xsonde&plev = plev
> xsonde!2 = "sonde"
>
> *printVarSummary*(xsonde)
>
>
> On Tue, Sep 3, 2019 at 3:51 PM Andrew Kren - NOAA Affiliate <
> andrew.kren at noaa.gov> wrote:
>
>> So what you are saying is that I should interpolate to the location
>> first. Then after I do that, I can interpolate multiple locations to the
>> time? This is how I understand it given the limitation.
>>
>> Thanks,
>>
>> On Tue, Sep 3, 2019 at 2:59 PM Dennis Shea <shea at ucar.edu> wrote:
>>
>>> I assume you are processing multiple GRIB2 files using *addfiles.*
>>> This is appropriate if the same grid is applicable to all files.
>>> *addfiles *will concatenate all the file times. You only need to use
>>> one grid.
>>> ---
>>> In the case of the HWRF grids (possibly) changing at each time step, you
>>> will have to step through the files manually. This means you will have to
>>> do some book-keeping.
>>>
>>> dirHWRF = "..."
>>> filHWRF = "..."
>>> pthHWRF = dirHWRF + filHWRF
>>> nfil = dimsizes(filHWRF)
>>>
>>> SLAT = ; one or more Sonde lat
>>> SLON =
>>>
>>> do nf=0,nfil-1
>>> f = addfile(pthHWRF(nf),"r")
>>> x = f->X
>>> xsonde(nf.....) = *linint2_points*(x&lon_0, x&lat_0,
>>> False,SLAT,SLON 0)
>>>
>>> end do
>>>
>>> Good Luck
>>>
>>> On Tue, Sep 3, 2019 at 1:19 PM Andrew Kren - NOAA Affiliate via ncl-talk
>>> <ncl-talk at ucar.edu> wrote:
>>>
>>>> Dear ncl-talk,
>>>>
>>>> I am working to interpolate a high-resolution model dataset (HWRF) to
>>>> the time and location of a dropsonde location. In the past for other
>>>> datasets, I have used linint1_n and linint2_points to interpolate in time
>>>> and space to the dropsonde position.
>>>>
>>>> However, HWRF, which is 6-hour analysis output, has a varying grid
>>>> every 6 hours as the domain "moves" with the tropical cyclone. Hence, the
>>>> lat/lon locations are not fixed in time. They have the same dimension sizes
>>>> but the points vary depending on the cyclone center. This led to an error
>>>> in my linint2_points function:
>>>>
>>>> fatal:linint2_points: If xi is not one-dimensional, then it must have
>>>> one less dimension than fi
>>>>
>>>>
>>>>
>>>> I know why this is happening. I read in two HWRF files around the
>>>> dropsonde time to do the interpolation. I used "join" to join the files and
>>>> my lat/lon arrays are (time, lat) and (time, lon). So I get the error.
>>>>
>>>>
>>>> But my question is how best to proceed? This is technically not an
>>>> unstructured grid, but could be considered one. Should I use ESMF?
>>>>
>>>>
>>>> Thanks,
>>>>
>>>> --
>>>> Andrew Kren
>>>> Assistant Scientist
>>>> University of Miami CIMAS - NOAA/AOML
>>>> Global Observing Systems Analysis (GOSA) Group
>>>> NOAA/AOML Quantitative Observing System Assessment Program (QOSAP)
>>>> _______________________________________________
>>>> ncl-talk mailing list
>>>> ncl-talk at ucar.edu
>>>> List instructions, subscriber options, unsubscribe:
>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>
>>>
>>
>> --
>> Andrew Kren
>> Assistant Scientist
>> University of Miami CIMAS - NOAA/AOML
>> Global Observing Systems Analysis (GOSA) Group
>> NOAA/AOML Quantitative Observing System Assessment Program (QOSAP)
>>
>
--
Andrew Kren
Assistant Scientist
University of Miami CIMAS - NOAA/AOML
Global Observing Systems Analysis (GOSA) Group
NOAA/AOML Quantitative Observing System Assessment Program (QOSAP)
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20190904/816837c9/attachment.html>
More information about the ncl-talk
mailing list