Re: re-converting 1D to 3D array.

Dave ... THX ... this *is* the way to do this.

D

On 03/10/2010 01:49 PM, Dave Allured wrote:
> Sho,
>
> My mistake, I did not notice the 2-D coordinates before. There is
> an efficient method to handle subsetting for 2-D coordinates. Use
> the region_ind function to get the bounding x and y subscripts for
> the region of interest:
>
> http://www.ncl.ucar.edu/Document/Functions/Contributed/region_ind.shtml
>
> Here is an untested upgrade to my first suggestion. This gets the
> bounding subscript values, then reads a 3-D subset data array of the
> desired dimensions.
>
> latS = 37
> latN = 47
> lonW = 269
> lonE = 279
>
> lat = fin->lat ; read 2-D coordinates
> lon = fin->lon
>
> ji = region_ind (lat, lon, latS, latN, lonW, lonE)
> j1 = ji(0) ; lat start ; get subset indices
> j2 = ji(1) ; lat last
> i1 = ji(2) ; lon start
> i2 = ji(3) ; lon last
>
> prc = fin->pr(:,{i1:i2},{j1:j2}) ; read 3-D subset
>
> lats_region = lat({i1:i2},{j1:j2}) ; make matching subset
> lons_region = lon({i1:i2},{j1:j2}) ; of 2-D coordinates
>
> printVarSummary (prc) ; check results before proceeding
> print (lats_region)
> print (lons_region)
> exit
>
> There may be be a few extra grid points along the edges where the
> subset data grid is out of alignment with a perfect rectangle as
> defined with lat/lon boundaries. You can then set these edge points
> to missing values if needed, but it is simpler for many applications
> to just leave them alone. I hope this is the right solution for you.
>
> --Dave A.
> NOAA/PSD/CIRES
>
> Sho Kawazoe wrote:
>> Hi,
>>
>> And thanks for the quick reply. However, based on your answer,
>> this assumes that the pr variable in the NetCDF file is originally in (time,
>> lat, lon) dimensions correct? It is in (time, xc, yc) coordinate, so is
>> there an maneuver to get around this?
>>
>> Sho
>>
>> On Tue, Mar 9, 2010 at 6:46 PM, Dave Allured<dave.allured@noaa.gov> wrote:
>>
>>> Sho Kawazoe,
>>>
>>> This is the second statement in your program:
>>>
>>> prc = fin->pr
>>>
>>> Move up the four lines for lat1, lat2, lon1, lon2; and try this instead:
>>>
>>> prc = fin->pr(:,{lon1:lon2},{lat1:lat2})
>>>
>>> printVarSummary(prc)
>>> print(prc&lon)
>>> print(prc&lat)
>>> exit
>>>
>>> If latitudes in the file go from high to low, then reverse the order
>>> of lat1 and lat2 above. Same for longitudes.
>>>
>>> Somehow you started making this reader more complicated than it
>>> should be. This should work better. We do not need to discuss
>>> conversion of dimensions unless your problem is more complicated
>>> than it looks.
>>>
>>> --Dave A.
>>> NOAA/PSD/CIRES
>>>
>>> Sho Kawazoe wrote:
>>>> NCL users,
>>>>
>>>> I'm still very raw in my knowledge with ncl, and was wondering if anyone
>>>> can help me out.
>>>>
>>>> I'm currently using precipitation data with MM5 (from EGS), and trying to
>>>> mask a specific lat and lon point, with the ultimate purpose of
>>> outputting
>>>> data within that specific grid, and truncating the data to where when I
>>> end
>>>> up concatenating other files, I am well below the 2gb threshold.
>>>>
>>>> Steps leading to retrieving the desiredvalues in 1D looks fine, and the
>>>> output file is well below the original input file size. However, when I
>>>> reconvert the 1D to 3D, the file is back to the same size as the input
>>> file,
>>>> and it looks like missing values are inputted into the output file, which
>>>> doesn't make the file any smaller.
>>>>
>>>> Is there a way to reconvert the 1D array to 3D array and keep only the
>>>> precip values we want, without the missing values(keep in mind that the
>>>> associated metadata needs to be there as well)? A simple onedtond doesn't
>>>> look like the solution, and I'm misunderstanding the coding.
>>>>
>>>> The script is as follows. Any help is appreciated.
>>>>
>>>> ;*****************************************************
>>>> ; Read original NetCDF file.
>>>> ;*****************************************************
>>>>
>>>> fin = addfile ("pr_MM5I_1979010103.nc","r")
>>>>
>>>> ;******************************************************
>>>> ; Upload netCDF file variables, set dimensions
>>>> ;******************************************************
>>>>
>>>> prc = fin->pr ; obtain pr(time,xc,yc) values. Dimension is
>>> in
>>>> 3D
>>>> lat = fin->lat
>>>> lon = fin->lon
>>>> time = fin->time
>>>> Lambert_Conformal = fin->Lambert_Conformal
>>>>
>>>> dim = dimsizes(prc)
>>>> ntime = dim(0)
>>>> nlat = dim(1)
>>>> nlon = dim(2)
>>>>
>>>> latlonmask2 = prc(:,:,:) ; assign dim size to masking var
>>>> latlonmask2 = 0 ; initialization
>>>>
>>>> ;**************************************************************
>>>> ; Indicate coordiate you want to analyze. Make sure you
>>>> ; check the NetCDF file for correct lat lon specification!!
>>>> ;
>>>> ; Mask through time, xc, and yc for entire precip variable to
>>>> ; analyze specifc regions
>>>> ;**************************************************************
>>>>
>>>> lat1 = 37
>>>> lat2 = 47
>>>> lon1 = 269
>>>> lon2 = 279
>>>>
>>>> do k = 1, ntime-1
>>>> do j = 1, nlat-1
>>>> do i = 1, nlon-1
>>>> if (lat(j,i).ge.lat1 .and. lat(j,i).lt.lat2 .and. lon(j,i).gt.lon1
>>>> .and. lon(j,i).le.lon2) then
>>>> latlonmask2(k,j,i) = 1
>>>> end if
>>>> end do
>>>> end do
>>>> end do
>>>>
>>>> desiredvalues = mask(prc, latlonmask2, 1)
>>>>
>>>> ;********************************************************************
>>>> ; Convert the masked region to a 1d array, which allows us to then
>>>> ; use the ncl function "ind" to keep only the data from the region
>>>> ; we are interested in looking at. This also keeps all the
>>>> ; metadata associated with each precip value.
>>>> ;*******************************************************************
>>>>
>>>> desired1d = ndtooned (desiredvalues) ; convert to 1d array
>>>> indices = ind (.not.ismissing(desired1d)) ; keep precip values that
>>>> have no missing data.
>>>> output1d = desired1d(indices) ; all non-missing values
>>>> desiredvalues = onedtond(output1d,dimsizes(desiredvalues))
>>>> delete (desired1d) ; clean up if you wish.
>>>> delete (indices)
>>>>
>>>> ;******************************************************************
>>>> ; Create a new nc output file
>>>> ;******************************************************************
>>>>
>>>> system("/bin/rm -f region_output.nc") ; remove any
>>> pre-existing
>>>> file
>>>> fout = addfile("region_output.nc","c")
>>>> filedimdef(fout,"time",-1,True) ; Makes time an UNLIMITED
>>>> dimension (as suggested by NCL)
>>>>
>>>> fout->prc = output1d
>>>>
>>>> fout->lat = lat
>>>> fout->lon = lon
>>>>
>>>> end
>>>>
>>>>
>>>> ------------------------------------------------------------------------
>>>>
>>>> _______________________________________________
>>>> 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
>>>
>>
>>
>>
>>
>> ------------------------------------------------------------------------
>>
>> _______________________________________________
>> 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

-- 
======================================================
Dennis J. Shea                  tel: 303-497-1361    |
P.O. Box 3000                   fax: 303-497-1333    |
Climate Analysis Section                             |
Climate & Global Dynamics Div.                       |
National Center for Atmospheric Research             |
Boulder, CO  80307                                   |
USA                        email: shea 'at' ucar.edu |
======================================================
_______________________________________________
ncl-talk mailing list
List instructions, subscriber options, unsubscribe:
http://mailman.ucar.edu/mailman/listinfo/ncl-talk