[ncl-talk] Filtering dates of a netcdf file (include number of gridpoints as additional criteria)

[Lyndz]

Dear Rashed,

Thank you so much for your help!
This works. I just need to tweak the threshold.

Sincerely,
Lyndz


On Tue, Sep 10, 2019 at 3:15 AM Rashed Mahmood <rashidcomsis at gmail.com>
wrote:

>
> I think your condition (b) would be the only condition that is needed,
> since if condition b is true then condition "a" must be true. So, it
> appears that you do not need condition "a".
> For example, if 80% of grid points for a given time and region has
> precipitation values greater than zero, then for that particular time and
> region area average value will always be above zero.
>
> So I modified your script (see attached) and added a loop to find out
> those time steps:
>
> ; condition b:
>
>    dims      = dimsizes(rain_region)
>    ntime     = dims(0)
>    tot_point = dims(1)*dims(2)*1.
>    ids_valid = new(ntime,integer)
>    req_perc  = 80. ; 40.
>
>    do n=0,ntime-1
>       tmp_rain = rain_region(n,:,:)
>       oneD     = ndtooned(tmp_rain)
>       npnts    = 1.*dimsizes(ind(oneD.gt.0))
>       perc     = 100.*(npnts/tot_point)
>
>       if(perc.gt.req_perc)then
>         ids_valid(n) = n
>        else
>         print(n+ " ... pertange of grid points is less than 80% = > "+perc)
>       end if
>       delete([/tmp_rain,oneD,npnts,perc/])
>    end do
>
>    if(.not.all(ismissing(ids_valid)))then
>        n_thres  = ids_valid(ind(.not.ismissing(ids_valid)))
>    else
>      print(" No valid time step for the specified condition: Exiting")
>     exit
>    end if
>
>
> You need to play with the variable    req_perc  = 80., by changing to
> lower values such 40. or anything else to see how it works. In your data
> file no time step would satisfy the 80% condition.
> hope that helps.
>
>
>
>
>
>
>
>
>
> On Mon, Sep 9, 2019 at 9:20 AM Lyndz <olagueralyndonmark429 at gmail.com>
> wrote:
>
>> Hi,
>>
>> This includes all time steps. But, in the file that I sent, I want to *count
>> the number of gridpoints in each time step* with rainfall above 0. I
>> have to include this in the time filter. Specifically, these two
>> conditions: the area averaged rainfall is > 0 *and* at least 80% of the
>> total gridpoints have rainfall > 0.
>>
>> Over the region I specified (rain_region), it has 40 x 30 x 10 (ntim x
>> nlat x nlon).
>> The above command you gave results to 5830.
>> I can get a similar result by just doing this: a = num(rain_region.gt.0)
>>
>> Sincerely,
>>
>>
>> On Tue, Sep 10, 2019 at 12:02 AM Rashed Mahmood <rashidcomsis at gmail.com>
>> wrote:
>>
>>> Yes, there can be several ways of counting grid points with rainfall
>>> above 0:
>>>
>>> e.g. say your variable name is "precp"
>>>
>>> oneD = ndtooned(precip)
>>>
>>>  nVal = dimsizes(ind(oneD.gt.0))
>>>
>>> nVal is the number of time precip is greater than 0. Does this makes
>>> sense, let me know if does not.
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>>
>>> On Mon, Sep 9, 2019 at 11:36 AM Lyndz via ncl-talk <ncl-talk at ucar.edu>
>>> wrote:
>>>
>>>> Dear NCL experts,
>>>>
>>>> I would like to filter the dates in a NetCDF file when the following
>>>> conditions are met:
>>>> (a) Area averaged rainfall is greater than 0.
>>>>
>>>> Area:
>>>>   latS = 12.5
>>>>   latN = 20.0
>>>>   lonL = 120.0
>>>>   lonR = 122.5
>>>>
>>>> (b) At least 80% of the grid points within the specified area satisfies
>>>> (a).
>>>>
>>>> Attached in this email is the script that I am currently working on.
>>>> The sample data can be accessed here:
>>>> https://www.dropbox.com/s/spfnwb4rhwjo90o/aphro_1979_AMJ.nc?dl=0
>>>>
>>>>
>>>> *[Problem]*
>>>> I was able to do the condition (a) but is there a way to count the
>>>> number of grid points with rainfall above 0 over the specified region? How
>>>> can I add this as another condition for the date filter?
>>>>
>>>> The time filter is like this:
>>>>
>>>>   rain_threshold = 0
>>>>   rain_aave=wgt_areaave_Wrap(rain_region,1.0,1.0,0)
>>>>   printVarSummary(rain_aave)
>>>>   n_thres  = ind(rain_aave .gt. rain_threshold) ;time indices
>>>>   print(n_thres)
>>>>   print("=====")
>>>>
>>>> I'll appreciate any guidance about this.
>>>>
>>>> Sincerely,
>>>> Lyndz
>>>>
>>>>
>>>>
>>>> _______________________________________________
>>>> ncl-talk mailing list
>>>> ncl-talk at ucar.edu
>>>> List instructions, subscriber options, unsubscribe:
>>>> http://mailman.ucar.edu/mailman/listinfo/ncl-talk
>>>
>>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mailman.ucar.edu/pipermail/ncl-talk/attachments/20190910/a01dbd8c/attachment.html>