[Thread Prev][Thread Next][Index]

speed problem now Re: [ferret_users] how to control domain of mode INTERP without awkward file writing/reading



Russ' samplexy approach worked (well for 20 out of 21 sites compared to file writing and reading), but it took over 2 minutes to do the sampling for 21 data sites, while the old awkward file write/read approach takes
1 or 2 seconds. So unless there's some other ideas, I'm giving up.

I do wonder why samplexy took so long.

I tried Ansley's variable loading approach but that didn't work. Maybe for a future upgrade, could ferret have a load function that forces evaluation of fields.

Cheers,

lev

Lev Tarasov -   Dept of Physics and Physical Oceanography,
		Memorial University of Newfoundland.
                email: lev@xxxxxx
                http://www.physics.mun.ca/~lev/
		Sabbatical Tel (519)-821-3555
		Fax (709)-864-8739

On Fri, 9 Aug 2013, Russ Fiedler wrote:


Hi Lev,

Ah, what you've done here is interpolate the (0,1) mask. Yes, this will create all sorts of problems as will interpolating quantities generated later with the mask.

I think what you need to do is sample/interpolate the values of H and Hw to the observed locations and then calculate the masks etc to do your calculation. So maybe something like

let h_sample=samplexy(h,obslon,obslat)
let hw_sample=samplexy(hw,obslon,obslat)

LET mask_grounded = if H_sample*rhoI/rhoO gt Hw_sample then 1 else 0

Cheers,
Russ




On 08/08/13 21:46, Lev Tarasov wrote:
Hi Russ;

I take it that H and Hw are on different grids (space?,time?) and that is why
you are using the set interpolate mode. I think the problem is that

H and Hw are on the same grid (see below). I also tried all the following 3 variants you suggested, and they had no impact:
 LET mask_grounded = if H*rhoI/rhoO gt Hw[g=H@AVE] then 1 else 0
 LET mask_grounded = if H*rhoI/rhoO gt Hw[g=H@LIN] then 1 else 0
 LET mask_grounded = if H*rhoI/rhoO gt Hw[g=H@MIN] then 1 else 0

I need set mode INTERP to interpolate from the 40km model grid to the location of the observational data. However, that causes the mask to
take fractional values which breaks the event crossing algorithm:

yes? stats mask_grounded
             IF H*RHOI/RHOO GT HW[G=H@MIN] THEN 1 ELSE 0
             LONGITUDE: 143.1E (interpolated)
             LATITUDE: 66.5S (interpolated)
             T: -120.05 to 0.05
 Total # of data points: 1201 (1*1*1*1201)
 Minimum value: 0
 Maximum value: 0.46
 Mean    value: 0.23324 (unweighted average)

So how can I keep set mode INTERP but have mask_grounded only
take on values 1 and 0?

Cheers,

Lev


yes? show grid Hw
    GRID (G001)
 name       axis              # pts   start                end
 XLONA1    LONGITUDE          360mr   0.5E                 0.5W
 YLATA5    LATITUDE            75 r   89.75S               52.75S
 normal    Z
 T120K100YR T                1201 r   -120                 0
yes? Show grid H
    GRID (G001)
 name       axis              # pts   start                end
 XLONA1    LONGITUDE          360mr   0.5E                 0.5W
 YLATA5    LATITUDE            75 r   89.75S               52.75S
 normal    Z
 T120K100YR T                1201 r   -120                 0


Lev Tarasov -   Dept of Physics and Physical Oceanography,
        Memorial University of Newfoundland.
                email: lev@xxxxxx
                http://www.physics.mun.ca/~lev/
        Sabbatical Tel (519)-821-3555
        Fax (709)-864-8739

On Thu, 8 Aug 2013, Russ Fiedler wrote:


Hi Lev,

I take it that H and Hw are on different grids (space?,time?) and that is why you are using the set interpolate mode. I think the problem is that interpolations are being performed by ferret in places and at times where you don't expect it.

I would explicitly regrid variables to a common grid and avoid any ambiguity.

Changing your mask definition to

LET mask_grounded = if H*rhoI/rhoO gt Hw[g=H@AVE] then 1 else 0
or
 LET mask_grounded = if H*rhoI/rhoO gt Hw[g=H@LIN] then 1 else 0
or maybe some other like @MIN may make sense.

That should allow you to run with interpolate mode off.


Cheers,
Russ



On 08/08/13 01:42, Lev Tarasov wrote:
I have a routine to figure out the time when a location last becomes free of grounded ice (ie for comparing model results to data for the last glacial cycle Antarctic ice sheet). To figure out whether the ice is floating, we need to use mode INTERP. However, this ends up messing up the @EVNT algorithm for finding out when the last crossing (in this case ungrounding of ice) occurs., ie: in the below, with mode interp on, mask_grounded (= if H*rhoI/rhoO gt Hw then 1 else 0) ends up with fractional values.

set region/x=143.14 /y=-66.48
! H and Hw are 2D space + 1D time fields
LET rhoI = 910  ![kg/m3]
LET rhoO = 1028 ![kg/m3]
LET mask_grounded = if H*rhoI/rhoO gt Hw then 1 else 0
LET grounded_event=mask_grounded[T=-120:0@EVNT:1]      ! find all the
instances when ice grounded (this starts at -120, incrementing for each event) LET last_event = grounded_event[T=@max] ! find the value of the last event LET mask_last = IF grounded_event EQ last_event THEN 1 ! make a mask based on last event value, ie domain with event value = last value LET idx= mask_last[T=@LOC:1] ! get the index of the first value in the mask \ie the last grounded event LET out = idx + 0.05 ! the 0.05 is to adjust the idx to output the time between grounded and floating

Now we can work around this by setting mode INTERP, saving the timeseries for each site, CAN DATA/ALL, can var/all, opening the saved timeseries and then doing the above masking operation. But for our large datasets, this
becomes awkward/slow. Anyway around this?

To be concrete, with the above work around for one data site and model, I get out= -31.75 without the work around I get out=-119.9 with CAN MODE INTERP and -70.57 with SET MODE INTERP

Cheers,

Lev

Lev Tarasov -   Dept of Physics and Physical Oceanography,
        Memorial University of Newfoundland.
                email: lev@xxxxxx
                http://www.physics.mun.ca/~lev/
        Sabbatical Tel (519)-821-3555
        Fax (709)-864-8739

This electronic communication is governed by the terms and conditions at
http://www.mun.ca/cc/policies/electronic_communications_disclaimer_2012.php




This electronic communication is governed by the terms and conditions at
http://www.mun.ca/cc/policies/electronic_communications_disclaimer_2012.php



This electronic communication is governed by the terms and conditions at
http://www.mun.ca/cc/policies/electronic_communications_disclaimer_2012.php


[Thread Prev][Thread Next][Index]
Contact Us
Dept of Commerce / NOAA / OAR / PMEL / Ferret

Privacy Policy | Disclaimer | Accessibility Statement