Trond, See if the attached file contains something close to what you need. It has no loops at all; I have not timed it, but it should be quite quick. I have given it only a cursory check, so I don't guarantee it works correctly. Depending on how your particular NetCDF interface works, you might need to copy arrays from NetCDF to ensure they are genuine ndarray objects. For plotting, you might want to try matplotlib. I think you will find it easier to use than GMT, especially if you are accustomed to matlab. Eric Trond Kristiansen wrote:
Hi again.
I have attached the function that the FOR loop is part of as a python file. What I am trying to do is to create a set of functions that will read the output files (NetCDF) from running the ROMS model (ocean model). The output file is organized in xi (x-direction), eta (y-direction), and s (z-direction) where the s-values are vertical layers and not depth. This particular function (z_slice) will find the closest upper and lower s-layer for a given depth in meters (e.g. -100). Then values from the two selcted layers will be interpolated to create a new layer at the selected depth (-100). The problem is that the s-layers follow the bathymetry and a particular s-layer will therefore sometimes be above and sometimes below the selected depth that we want to interpolate to. That's why I need a quick script that searches all of the layers and find the upper and lower layers for a given depth value (which is negative). The z_r is a 3D array (s,eta,xi) that is created using the netcdf module.
The main goal of these set of functions is to move away from using matlab, but also to speed things up. The sliced data array will be plotted using GMT or pyNGL.
Thanks for helping me. Cheers, Trond