Difference between revisions of "Matlab Processing Scripts"

From WikiROMS
Jump to navigationJump to search
Line 76: Line 76:


==Contact Regions and Contact Points Scripts==
==Contact Regions and Contact Points Scripts==
<wikitex>
These functions process the contact points for each contact region. It is an important step in configuring ROMS nesting applications. The objective is to create the contact points NetCDF file to be read by ROMS.  These functions are generic and have been tested extensively.  However, we recommend users to analyzed the data in the NetCDF file and activate the plotting switch to examine the contact point.  The plot maybe very dense when the connectivity involves a lot of files.  We need to zoom the plotting frames until there are readable and think carefully.  A expertise in ROMS is required to analyze the plots.
These function use the excellent Matlab function <span class="blue">inpolygon</span> to find if a point (X,Y) is inside a grid perimeter, on the grid perimeter edges, or outside the grid perimeter. In Matlab, we need to account for the zero shift of ROMS indices.  The logic is tricky and not that intuitive.  We need to take into account ROMS (I,J) indices and fractional coordinates ($\xi$, $\eta$).  For example, in Matlab the Imin and Jmin indices correspond to fractional coordinates:
$$\xi_{\psi}(I_{\psi}) = \xi_{\rho}(I_{\rho})$$
$$\eta_{\psi}(J_{\psi}) = \eta_{\rho}(J_{\rho})$$
since $I_{\psi} = I_{\rho}$ and $J_{\psi} = J_{\rho}$.  The number of points in each direction is $L_{\rho} = L_{\psi}+1$ and $M_{\rho} = M_{\psi}+1$, respectively.  Therefore, the zero shift of &rho;-points indices is done at the end.
The main scripts to carry this tasks are located in the '''matlab/grid''' sub-directory:
;<span id="c_contact"></span><span class="blue">c_contact.m</span>
:Creates a new contact coints NetCDF file for a ROMS nesting application. The contact point data is written elsewhere using [[Matlab_Processing_Scripts#write_contact|write_contact.m]]. if the <span class="green">Ndatum</span> argument is not provided, the '''datum''' dimension is is '''unlimited'''.<br />
::<span class="red">c_contact</span> &nbsp;(<span class="green">ncname</span>, <span class="green">spherical</span>, <span class="green">Ngrids</span>, <span class="green">Ndatum</span>)
:'''On Input:'''
::<span class="green">ncname</span>: &nbsp;  Ouptut contact points NetCDF file name (string)
::<span class="green">spherical</span>: &nbsp; Spherical grid switch (logical)
::<span class="green">Ngrids</span>: &nbsp; Number of nested grids
::<span class="green">Ndatum</span>: &nbsp; Total number of contact points (optional)
;<span id="contact"></span><span class="blue">contact.m</span>
:Sets ROMS nested grid contact points for each contact region. The order of the nested grid file names in the input cell array (Gnames) is important.  Set the file names in the order of nesting layers and time-stepping in ROMS. The <span class="green">MaskInterp</span> switch is only relevant when processing the contact region of a refinement grid. The setting of land/sea masking in the contact region is critical. Usually it is better to linearly interpolate the finer grid &psi;- u-, and v-masks from coarser grid than recomputing via [[uvp_masks]]. The recomputed masks may be one point smaller or larger. Use the <span class="green">MaskInterp</span> switch to either interpolate (true) or use [[uvp_masks]] (false).  Recall that we are not modifying the original refined grid mask, but just computing the mask in the contact region adjacent to finer grid from the coarser grid mask. This is only relevant when there is a land/sea masking features in any of the refinement grid physical boundary. If so, the user just need to experiment with <span class="green">MaskInterp</span> and edit such points during post-processing.<br />
::[<span class="green">S</span>, <span class="green">G</span>] = <span class="red">contact</span> &nbsp;(<span class="green">Gnames</span>, <span class="green">Cname</span>, <span class="green">Lmask</span>, <span class="green">MaskInterp</span>, <span class="green">Lplot</span>)
:'''On Input:'''
::<span class="green">Gnames</span>: &nbsp; Input Grid NetCDF full path file names  (cell array)
::<span class="green">Cname</span>: &nbsp; Ouptut contact points NetCDF file name (string)
::<span class="green">Lmask</span>: &nbsp; Switch to remove contact points over land (logical, optional, default false)
::<span class="green">MaskInterp</span>: &nbsp; Switch to interpolate &psi;-, u- and v-masks (true) or compute from interpolated &rho;-mask (false). The default is false (logical, optional).
::<span class="green">Cname</span>: &nbsp; Switch to plot various contact points figures (logical, optional, default true)
:'''On Ouput:'''
::<span class="green">S</span>: &nbsp; Nested grids contact points structure (struct array)
::<span class="green">G</span>: &nbsp; Nested grids structure (1 x Ngrids struct array)
;<span id="disp_contact"></span><span class="blue">disp_contact.m</span>
:Prints the nested grids contact points unique values for each Contact Region. The contact points structure is rich and this simple function facilitates to summarize the contact points values.<br />
::<span class="red">disp_contact</span> &nbsp;(<span class="green">S</span>)
:'''On Input:'''
::<span class="green">ncname</span>: &nbsp;  Contact points structure(struct array)
</wikitex>


==Initial Condition Scripts==
==Initial Condition Scripts==


==Support Scripts==
==Support Scripts==

Revision as of 01:22, 21 April 2012

Matlab Processing Scripts for Nesting

The Matlab processing scripts for nesting applications are divided in four categories:

  • Scripts to aid in the building of nested grids and associated NetCDF files,
  • Scripts to process contact regions and contact points and build its associated NetCDF file,
  • Scripts to process initial conditions data and NetCDF files, and
  • Generic support scripts.

Note Notice: Starting Matlab version 2012a, released on Feb 9, 2012, the native interface to NetCDF is the preferred method for processing NetCDF data in the scripts distributed in the ROMS repository matlab: https://www.myroms.org/svn/src/matlab. The native interface was introduced in Matlab version 2008b for NetCDF-3 type files. The NetCDF-4 support started in version 2010b. The support for HDF5 files was completed in version 2011a. The OpenDAP support began in version 2012a. If your Matlab version is older than 2008b, we highly recommend to update to the newest version. However, in the basic generic scripts we have switches for older versions to activate either the MEXNC interface for standard NetCDF files and the SNCTOOLS interface to process NetCDF files on an OpenDAP server.

Warning Caution: The ROMS Grid structure used in all the scripts described here is the one returned from by the script get_roms_grid.m.

Nested Grid Generation Scripts

Currently, it is not our intention to develop a fancy GUI tool to prepare ROMS nested grids. The capability of ROMS nesting is unique and requires special tools since we support mosaics, composite, and refinement grids in a variety of nested classes. Perhaps in the future an expert user will help us develop this GUI. However, the scripts described here are basic and fundamental. They can be used in the future to build a fancier interface.

Two strategies are possible when configuring nested grids for a particular application:

  1. Build an intermediate large coarse resolution grid and extract the composite and/or refinement grids from it using the functions dicussed below. The disadvantage of this approach is that we need work independently each grid for adjustments in the land/sea masking, target smoothed bathymetry, and matching domain volumes. Additionally, the finer resolution grids bathymetry can be under- or over-estimated due to the original coarser bathymetry. There is a linear definition of topographic features that are unresolved in the finer grids.
  2. Build an intermediate large very fine resolution grid, work on the desired land/sea masking, and target smoothed bathymetry. The intermediate grid must resolve the smaller resolution desired for a particular application. Then, use the scripts described below to extract the composite and/or refinement grids. We think that this strategy is advantageous in all aspects. However, it is cumbersome because we are dealing with a large grid. The land/sea masking maybe tedious. Anyway, the function landsea.m can be used to automatize the land/sea masking , but it still needs some hand editing using editmask.m.

The main scripts to carry this tasks are located in the matlab/grid sub-directory:

c_grid.m
Creates a new ROMS grid NetCDF file or modifies a existing ROMS grid NetCDF file.
status = c_grid  (Lp, Mp, Gname, NewFile, spherical)
On Input:
Lp:   Number of ρ-point in the ξ-direction
Mp:   Number of ρ-point in the η-direction
Gname:   Grid NetCDF file name (string)
NewFile:   Switch to create a new file (logical, optional, default false)
spherical:   Spherical grid switch (logical, optional, default true)
On Ouput:
status:   Error flag
coarse2fine
Given a coarse resolution ROMS Grid NetCDF file (Ginp), this function creates a finer resolution grid in the region specified by the coarser grid coordinates (Imin, Jmin) and (Imax, Jmax). Notice that (Imin, Jmin), and (Imax, Jmax) indices are in terms of the ψ-points because it actually define the physical boundaries of the refined grid. The grid refinement coefficient is specified with Gfactor.
F = coarse2fine  (Ginp, Gout, Gfactor, Imin, Imax, Jmin, Jmax)
On Input:
Ginp:   Input coarser Grid NetCDF file name (string)
Gout:   Output finer Grid NetCDF file name (string)
Gfactor:   Grid refinement factor (3,5,7,9,11,13,15,...)
Imin:   Coarse grid lower-left I-coordinate (ψ-point)
Imax:   Coarse grid upper-right I-coordinate (ψ-point)
Jmin:   Coarse grid lower-left J-coordinate (ψ-point)
Jmax:   Coarse grid upper-right J-coordinate (ψ-point)
On Ouput:
F:   Fine resolution grid structure (struct array)
fine2coarse
Given a finer resolution ROMS Grid NetCDF file (Ginp), this function creates a coarser resolution grid in the region specified by the finer grid coordinates (Imin, Jmin) and (Imax, Jmax). Notice that (Imin, Jmin), and (Imax, Jmax) indices are in terms of the ψ-points because it actually define the physical boundaries of the coarser grid. The grid coarseness coefficient is specified with Gfactor.
C = fine2coarse  (Ginp, Gout, Gfactor, Imin, Imax, Jmin, Jmax)
On Input:
Ginp:   Input finer Grid NetCDF file name (string)
Gout:   Output coarser Grid NetCDF file name (string)
Gfactor:   Grid coarseness factor (3,5,7,9,11,13,15,...)
Imin:   Finer grid lower-left I-coordinate (ψ-point)
Imax:   Finer grid upper-right I-coordinate (ψ-point)
Jmin:   Finer grid lower-left J-coordinate (ψ-point)
Jmax:   Finer grid upper-right J-coordinate (ψ-point)
On Ouput:
C:   Coarser resolution grid structure (struct array)
grid_extract
Given a larger ROMS Grid NetCDF file (Ginp), this function extracts and creates a sub-domain Grid NetCDF file (Gout) for the sampled region in terms of the larger grid coordinates (Imin, Jmin) and (Imax, Jmax). Notice that the (Imin, Jmin), and (Imax, Jmax) indices are located at ψ-points. They actually define the physical boundaries of the smaller sub-domain grid.
S = grid_extract  (Ginp, Gout, Imin, Imax, Jmin, Jmax)
On Input:
Ginp:   Input coaser Grid NetCDF file name (string) or an existing grid structure (struct array)
Gout:   Output smaller Grid NetCDF file name (string)
Imin:   Larger grid lower-left I-coordinate (ψ-point)
Imax:   Larger grid upper-right I-coordinate (ψ-point)
Jmin:   Larger grid lower-left J-coordinate (ψ-point)
Jmax:   Larger grid upper-right J-coordinate (ψ-point)
On Ouput:
S:   Smaller grid structure (struct array)

Contact Regions and Contact Points Scripts

<wikitex> These functions process the contact points for each contact region. It is an important step in configuring ROMS nesting applications. The objective is to create the contact points NetCDF file to be read by ROMS. These functions are generic and have been tested extensively. However, we recommend users to analyzed the data in the NetCDF file and activate the plotting switch to examine the contact point. The plot maybe very dense when the connectivity involves a lot of files. We need to zoom the plotting frames until there are readable and think carefully. A expertise in ROMS is required to analyze the plots.

These function use the excellent Matlab function inpolygon to find if a point (X,Y) is inside a grid perimeter, on the grid perimeter edges, or outside the grid perimeter. In Matlab, we need to account for the zero shift of ROMS indices. The logic is tricky and not that intuitive. We need to take into account ROMS (I,J) indices and fractional coordinates ($\xi$, $\eta$). For example, in Matlab the Imin and Jmin indices correspond to fractional coordinates:

$$\xi_{\psi}(I_{\psi}) = \xi_{\rho}(I_{\rho})$$ $$\eta_{\psi}(J_{\psi}) = \eta_{\rho}(J_{\rho})$$

since $I_{\psi} = I_{\rho}$ and $J_{\psi} = J_{\rho}$. The number of points in each direction is $L_{\rho} = L_{\psi}+1$ and $M_{\rho} = M_{\psi}+1$, respectively. Therefore, the zero shift of ρ-points indices is done at the end.

The main scripts to carry this tasks are located in the matlab/grid sub-directory:

c_contact.m
Creates a new contact coints NetCDF file for a ROMS nesting application. The contact point data is written elsewhere using write_contact.m. if the Ndatum argument is not provided, the datum dimension is is unlimited.
c_contact  (ncname, spherical, Ngrids, Ndatum)
On Input:
ncname:   Ouptut contact points NetCDF file name (string)
spherical:   Spherical grid switch (logical)
Ngrids:   Number of nested grids
Ndatum:   Total number of contact points (optional)
contact.m
Sets ROMS nested grid contact points for each contact region. The order of the nested grid file names in the input cell array (Gnames) is important. Set the file names in the order of nesting layers and time-stepping in ROMS. The MaskInterp switch is only relevant when processing the contact region of a refinement grid. The setting of land/sea masking in the contact region is critical. Usually it is better to linearly interpolate the finer grid ψ- u-, and v-masks from coarser grid than recomputing via uvp_masks. The recomputed masks may be one point smaller or larger. Use the MaskInterp switch to either interpolate (true) or use uvp_masks (false). Recall that we are not modifying the original refined grid mask, but just computing the mask in the contact region adjacent to finer grid from the coarser grid mask. This is only relevant when there is a land/sea masking features in any of the refinement grid physical boundary. If so, the user just need to experiment with MaskInterp and edit such points during post-processing.
[S, G] = contact  (Gnames, Cname, Lmask, MaskInterp, Lplot)
On Input:
Gnames:   Input Grid NetCDF full path file names (cell array)
Cname:   Ouptut contact points NetCDF file name (string)
Lmask:   Switch to remove contact points over land (logical, optional, default false)
MaskInterp:   Switch to interpolate ψ-, u- and v-masks (true) or compute from interpolated ρ-mask (false). The default is false (logical, optional).
Cname:   Switch to plot various contact points figures (logical, optional, default true)
On Ouput:
S:   Nested grids contact points structure (struct array)
G:   Nested grids structure (1 x Ngrids struct array)
disp_contact.m
Prints the nested grids contact points unique values for each Contact Region. The contact points structure is rich and this simple function facilitates to summarize the contact points values.
disp_contact  (S)
On Input:
ncname:   Contact points structure(struct array)


</wikitex>

Initial Condition Scripts

Support Scripts