Ocean Modeling Discussion

ROMS/TOMS

Search for:
It is currently Sun Jul 21, 2019 9:22 pm




Post new topic Reply to topic  [ 32 posts ] 

All times are UTC

Author Message
PostPosted: Mon Jun 17, 2019 5:08 pm 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
Hi everyone, please help me to create boundary and initial file from seagrid and GLBa0.08 - Hycom data.

I created seagrid nc file for my domain and have Hycom data. I also install pyroms_toolbox. My problem is I have no idea about how to make boundary and initial file. Could anybody help me a example by using pyroms_toolbox to do that please? especially how to fill in the file to run

My questions are

1) in seagrid how can I chose the number of vertical cell? (if I want my grid size is x=30; y=30; z=15)

2) for Grid_HYCOM, get_nc_Grid_HYCOM and make_remap_grid_file which information I should use in?
------------------------------------------------------------------------------------
class Grid_HYCOM(object):
"""
Grid object for HYCOM
"""

def __init__(self, lon_t, lat_t, lon_vert, lat_vert, mask_t, z_t, h, angle, name):

self.name = name

self.lon_t = lon_t
self.lat_t = lat_t

self.lon_vert = lon_vert
self.lat_vert = lat_vert

self.mask_t = mask_t

self.z_t = z_t

self.h = h

self.angle = angle
-----------------------------------------------------------------------------
3) should I run only one of flood, flood_fast and flood_fast_weighted or all of them?

Regard.


Top
 Profile  
Reply with quote  
PostPosted: Mon Jun 17, 2019 6:25 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
I have successfully made ROMS IC and BC files from HYCOM using the scripts in pyroms/examples/Arctic_HYCOM/. First you run make_remap_weights_file, then make_ic_file/make_bdry_file.
A more equatorial example is in Palau_HYCOM.

For the vertical, set the number of vertical points in the gridid.txt file. My ARCTIC4 chunk looks like:
Code:
id      = ARCTIC4
name    = ARCTIC4
grdfile = /import/AKWATERS/kshedstrom/gridpak/Arctic2/grid_Arctic_4.nc
N       = 50
grdtype = roms
Vtrans  = 4
theta_s = 7
theta_b = 2
Tcline  = 250


Top
 Profile  
Reply with quote  
PostPosted: Tue Jun 18, 2019 4:48 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
. First you run make_remap_weights_file, then make_ic_file/make_bdry_file.
A more equatorial example is in Palau_HYCOM.


Dear Kate

I download Hycom data with empty z (to import all z to the netcdf file)
then try to run make_remap_weights_file and got problem with this
------------------------------------------------------------------
Code:
>>> srcgrd = pyroms_toolbox.Grid_HYCOM.get_nc_Grid_HYCOM('/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/Hycomdata/1906201200.nc')

Traceback (most recent call last):
  File "<pyshell#6>", line 1, in <module>
    srcgrd = pyroms_toolbox.Grid_HYCOM.get_nc_Grid_HYCOM('/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/Hycomdata/1906201200.nc')
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/Grid_HYCOM/get_nc_Grid_HYCOM.py", line 19, in get_nc_Grid_HYCOM
    depth = nc.variables['z'][:]
KeyError: 'z'

------------------------------------------------------------------
by open the hycom netcdf file data by IDV and I have all z data can be visualized in 3d model so I do not understand why nc.variables['z'][:] have that error?

by running ncdump -h /work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/Hycomdata/1906201200.nc > my_file.cdl I got 40 depth levels like this.
-----------------------------------------------------------------
Code:
netcdf \1906201200 {
dimensions:
   time = 5 ;
   lat = 626 ;
   lon = 313 ;
   depth = 40 ;
variables:
...
float salinity(time, depth, lat, lon) ;
      salinity:_CoordinateAxes = "time_run time depth lat lon" ;
      salinity:units = "psu" ;
      salinity:long_name = "Salinity" ;
      salinity:standard_name = "sea_water_salinity" ;
      salinity:NAVO_code = 16 ;
      salinity:coordinates = "time_run time depth lat lon " ;
   double depth(depth) ;
      depth:units = "m" ;
      depth:long_name = "Depth" ;
      depth:standard_name = "depth" ;
      depth:positive = "down" ;
      depth:axis = "Z" ;
      depth:NAVO_code = 5 ;
      depth:_CoordinateAxisType = "Height" ;
      depth:_CoordinateZisPositive = "down" ;
...

----------------------------------------------------------------
my file /get_nc_Grid_HYCOM.py", line 19 got
Code:
    depth = nc.variables['Z'][:]
#   depth = nc.variables['z'][:]


I tried booth Z and z with same error
-----------------------------------------------------------------
Thank you so much for your allway very quick and helpful replies.

Regard.


Top
 Profile  
Reply with quote  
PostPosted: Tue Jun 18, 2019 5:20 am 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Looks like it should be 'depth', not 'z'.


Top
 Profile  
Reply with quote  
PostPosted: Tue Jun 18, 2019 6:05 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
Looks like it should be 'depth', not 'z'.


Dear Kate
you are right I changed it to depth and pass, and have this error after that

Code:
>>> srcgrd = pyroms_toolbox.Grid_HYCOM.get_nc_Grid_HYCOM('/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/Hycomdata/1906201200.nc')

Traceback (most recent call last):
  File "<pyshell#4>", line 1, in <module>
    srcgrd = pyroms_toolbox.Grid_HYCOM.get_nc_Grid_HYCOM('/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/Hycomdata/1906201200.nc')
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/Grid_HYCOM/get_nc_Grid_HYCOM.py", line 21, in get_nc_Grid_HYCOM
    var = nc.variables['sea_water_temperature'][0,:,1:-1,1:-1]
KeyError: 'sea_water_temperature'


my cdl file have this information

Code:
float water_temp(time, depth, lat, lon) ;
      water_temp:_CoordinateAxes = "time_run time depth lat lon" ;
      water_temp:units = "degC" ;
      water_temp:long_name = "Water Temperature" ;
      water_temp:standard_name = "sea_water_temperature" ;
      water_temp:NAVO_code = 15 ;
      water_temp:comment = "in-situ temperature" ;
      water_temp:coordinates = "time_run time depth lat lon " ;


So what should I need to replace for "sea_water_temperature" ?
Change it to water_temp and passed

then I got this
Code:
>>> srcgrd = pyroms_toolbox.Grid_HYCOM.get_nc_Grid_HYCOM('/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/Hycomdata/1906201200.nc')

Traceback (most recent call last):
  File "<pyshell#4>", line 1, in <module>
    srcgrd = pyroms_toolbox.Grid_HYCOM.get_nc_Grid_HYCOM('/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/Hycomdata/1906201200.nc')
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/Grid_HYCOM/get_nc_Grid_HYCOM.py", line 24, in get_nc_Grid_HYCOM
    lon_t = lon[1:-1,1:-1]
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/numpy/ma/core.py", line 3174, in __getitem__
    dout = self.data[indx]
IndexError: too many indices for array


My installed pyroms is for python 2.7 with included Hycom folder with scripts so I hope these scripts are wrote for python 2.7 but it seem that "IndexError: too many indices for array" related to version conflict doesn't it?


Thank you
Regard.


Top
 Profile  
Reply with quote  
PostPosted: Tue Jun 18, 2019 7:13 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
No, I don't think you have a version problem, unless it's a HYCOM version problem.
Code:
lon_t = lon[1:-1,1:-1]
Longitude is often a 1-D array, not a 2-D array. What do you have?


Top
 Profile  
Reply with quote  
PostPosted: Wed Jun 19, 2019 2:42 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
No, I don't think you have a version problem, unless it's a HYCOM version problem.
Code:
lon_t = lon[1:-1,1:-1]
Longitude is often a 1-D array, not a 2-D array. What do you have?


Dear Kate
I am sorry, I still do not understand my problem.
I would like to make IC BC from this data https://ncss.hycom.org/thredds/ncss/grid/GLBy0.08/latest/dataset.html
and I have get_nc_Grid_HYCOM.py code after install pyroms like this

Code:
import numpy as np
import pyroms
import netCDF4
from mpl_toolkits.basemap import pyproj
from pyroms_toolbox.Grid_HYCOM import Grid_HYCOM


def get_nc_Grid_HYCOM(grdfile, name='GLBa0.08_NEP'):

    """
    grd = get_nc_Grid_HYCOM(grdfile)

    Load grid object for HYCOM_GLBa0.08_NEP
    """

    nc = netCDF4.Dataset(grdfile)
    lon = nc.variables['lon'][:]
    lat = nc.variables['lat'][:]
    depth = nc.variables['z'][:]
    var = nc.variables['temp'][0,:,1:-1,1:-1]
    nc.close()

    lon_t = lon[1:-1,1:-1]
    lat_t = lat[1:-1,1:-1]

    lon_vert = 0.5 * (lon[:,1:] + lon[:,:-1])
    lon_vert = 0.5 * (lon_vert[1:,:] + lon_vert[:-1,:])

    lat_vert = 0.5 * (lat[1:,:] + lat[:-1,:])
    lat_vert = 0.5 * (lat_vert[:,1:] + lat_vert[:,:-1])

    mask_t = np.array(~var[:].mask, dtype='int')

    z_t = np.tile(depth,(mask_t.shape[2],mask_t.shape[1],1)).T

    depth_bnds = np.zeros(len(depth)+1)
    for i in range(1,len(depth)):
        depth_bnds[i] = 0.5 * (depth[i-1] + depth[i])
    depth_bnds[-1] = 5750

    bottom = pyroms.utility.get_bottom(var[::-1,:,:], mask_t[0], spval=var.fill_value)
    nlev = len(depth)
    bottom = (nlev-1) - bottom
    h = np.zeros(mask_t[0,:].shape)
    for i in range(mask_t[0,:].shape[1]):
        for j in range(mask_t[0,:].shape[0]):
            if mask_t[0,j,i] == 1:
                h[j,i] = depth_bnds[int(bottom[j,i])+1]


    geod = pyproj.Geod(ellps='WGS84')
    az_forward, az_back, dx = geod.inv(lon_vert[:,:-1], lat_vert[:,:-1], lon_vert[:,1:], lat_vert[:,1:])
    angle = 0.5 * (az_forward[1:,:] + az_forward[:-1,:])
    angle = (90 - angle) * np.pi/180.


    return Grid_HYCOM(lon_t, lat_t, lon_vert, lat_vert, mask_t, z_t, h, angle, name)



Thank you
Regard


Top
 Profile  
Reply with quote  
PostPosted: Wed Jun 19, 2019 5:10 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Are you sure you are downloading on the GLBa0.08 HYCOM grid? They have moved that to the "inactive" status. It was on a tripole grid and needed lat,lon to be functions if both i and j. The newer HYCOM products are on a simpler grid, where lat is a function of j only and lon is a function of i only. I'm afraid I haven't tried to run these scripts with the new HYCOM products. It's on my list of things to do, just not necessarily this week.


Top
 Profile  
Reply with quote  
PostPosted: Sat Jun 22, 2019 2:26 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
Are you sure you are downloading on the GLBa0.08 HYCOM grid? They have moved that to the "inactive" status. It was on a tripole grid and needed lat,lon to be functions if both i and j. The newer HYCOM products are on a simpler grid, where lat is a function of j only and lon is a function of i only. I'm afraid I haven't tried to run these scripts with the new HYCOM products. It's on my list of things to do, just not necessarily this week.


Dear Kate
Thank you for your suggestion. I double check and found out I have downloaded wrong data for GLBa0.08. Actually I would like to have that data for predict/ forecasting works but at the moment, for just practice using ROMS, It is better for me if I keep up with GLBa0.08 examples to understand how pyroms work.

I followed instruction in this topic https://www.myroms.org/forum/viewtopic.php?f=23&t=4819
mdessert wrote
Quote:
1-I've well downloaded files from HYCOM.org.
2-I've formatted them through get_hycom_GLBa0.08_'varname'_2014.py.
3-I've built my source grid through get_hycom_GLBa0.08_grid.py.
4-I've computed my weigth files from make_remap_weights_file.py.
5-And then, I've 'concatenated' all my variable files to obtain one big file HYCOM_GLBa0.08_2014_095.nc.

at 1- I am not sure how he downloaded these files, I try to run 2- then get the *.nc file for salt, ssh, temp, u and v

then I run 3- once time and get the GLBa0.08 grid.nc file
when I run 4- for the 1st time I stuck as mention bellow.
I do not know how to run step 5- to have all in one big file HYCOM_GLBa0.08_2014_095.nc


I stuck at running make_remap_weights_file.py :

Code:
Traceback (most recent call last):
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/make_remap_weights_file.py", line 7, in <module>
    srcgrd = pyroms_toolbox.Grid_HYCOM.get_nc_Grid_HYCOM('/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/HYCOM_GLBa0.08_PALAU_grid.nc')
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/Grid_HYCOM/get_nc_Grid_HYCOM.py", line 52, in get_nc_Grid_HYCOM
    h[j,i] = depth_bnds[bottom[j,i]+1]
IndexError: only integers, slices (`:`), ellipsis (`...`), numpy.newaxis (`None`) and integer or boolean arrays are valid indices


I am sorry because I am not a code/programmer then It should take me time to pass basic questions like this.
Regard and have a nice weekend.


Top
 Profile  
Reply with quote  
PostPosted: Sat Jun 22, 2019 3:35 am 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
My latest version of that troublesome line is:
Code:
h[j,i] = depth_bnds[int(bottom[j,i])+1]
Sometimes it pays to "git pull".

Step 5 happens with this script:
Code:
#!/bin/ksh

print -n "Enter year to process: "; read year; print ""

leap=`echo $(($year % 4))`

if [ $leap == 0 ] ; then
    nday=366
else
    nday=365
fi
#nday=31

#set -A days {1..331}
#set -A days {2..4}
#set -A days {33..$nday}
set -A days {1..$nday}

for day in ${days[@]} ; do
    day=`echo $day | awk '{printf "%03d", $1}'`
    ncks -a -O HYCOM_GLBa0.08_ssh_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    ncks -a -A HYCOM_GLBa0.08_temp_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    ncks -a -A HYCOM_GLBa0.08_salt_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    ncks -a -A HYCOM_GLBa0.08_u_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    ncks -a -A HYCOM_GLBa0.08_v_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    rm -f HYCOM_GLBa0.08_ssh_${year}_${day}.nc HYCOM_GLBa0.08_temp_${year}_${day}.nc \ HYCOM_GLBa0.08_salt_${year}_${day}.nc HYCOM_GLBa0.08_u_${year}_${day}.nc \
HYCOM_GLBa0.08_v_${year}_${day}.nc
done
Gosh, should convert to bash.


Top
 Profile  
Reply with quote  
PostPosted: Sat Jun 22, 2019 3:56 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
My latest version of that troublesome line is:
Code:
h[j,i] = depth_bnds[int(bottom[j,i])+1]
Sometimes it pays to "git pull".



Dear Kate
by looking in get_nc_Grid_HYCOM.py I have

Code:
for i in range(1,len(depth)):
        depth_bnds[i] = 0.5 * (depth[i-1] + depth[i])
    depth_bnds[-1] = 5750


It suggested that after that, depth_bnds become not integer anymore (0.5*....), then could it be leading to the "IndexError: only integers, slices (`:`), ellipsis (`...`), numpy.newaxis (`None`) and integer or boolean arrays are valid indices" error?
Regard.


Top
 Profile  
Reply with quote  
PostPosted: Sat Jun 22, 2019 4:03 am 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Sometimes Python evolves so that something that used to run no longer runs. Sometimes when that happens I can figure out how to change it and get it pushed to github. You need to update your code through either "git pull" or just editing it in place.


Top
 Profile  
Reply with quote  
PostPosted: Sat Jun 22, 2019 8:55 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
Sometimes Python evolves so that something that used to run no longer runs. Sometimes when that happens I can figure out how to change it and get it pushed to github. You need to update your code through either "git pull" or just editing it in place.


Dear Kate
Thank you so much, I do git pull and updated the code, then pyroms go to this

Code:
Traceback (most recent call last):
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/make_remap_weights_file.py", line 8, in <module>
    dstgrd = pyroms.grid.get_ROMS_grid('PALAU1')
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms/grid.py", line 467, in get_ROMS_grid
    gridinfo = ROMS_gridinfo(gridid,hist_file=hist_file,grid_file=grid_file)
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms/grid.py", line 71, in __init__
    self._get_grid_info(grid_file,hist_file)
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms/grid.py", line 82, in _get_grid_info
    data = open(gridid_file,'r')
TypeError: coercing to Unicode: need string or buffer, NoneType found


I think it is related to gridid.txt. My gridid.txt has no infor about PALAU1 (I used script in PALAU_HYCOM in example folder. So how can I find information to fill in gridid.txt to run?

Regard.


Top
 Profile  
Reply with quote  
PostPosted: Sat Jun 22, 2019 3:40 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
I can give you the PALAU1 chunk of gridid.txt, but then it points to a grid file which you don't have either. Why don't you put in a chunk for your grid, since that is what you want anyway?
Code:
id      = PALAU1
name    = PALAU1
grdfile = /import/AKWATERS/kshedstrom/gridpak/Palau/grid_Palau_1.nc
N       = 50
grdtype = roms
Vtrans  = 4
theta_s = 7
theta_b = 2
Tcline  = 250
This is the PALAU chunk which has the vertical coordinates we use for most domains these days.


Top
 Profile  
Reply with quote  
PostPosted: Sun Jun 23, 2019 3:34 pm 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
My latest version of that troublesome line is:
Code:
h[j,i] = depth_bnds[int(bottom[j,i])+1]
Sometimes it pays to "git pull".

Step 5 happens with this script:
Code:
#!/bin/ksh

print -n "Enter year to process: "; read year; print ""

leap=`echo $(($year % 4))`

if [ $leap == 0 ] ; then
    nday=366
else
    nday=365
fi
#nday=31

#set -A days {1..331}
#set -A days {2..4}
#set -A days {33..$nday}
set -A days {1..$nday}

for day in ${days[@]} ; do
    day=`echo $day | awk '{printf "%03d", $1}'`
    ncks -a -O HYCOM_GLBa0.08_ssh_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    ncks -a -A HYCOM_GLBa0.08_temp_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    ncks -a -A HYCOM_GLBa0.08_salt_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    ncks -a -A HYCOM_GLBa0.08_u_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    ncks -a -A HYCOM_GLBa0.08_v_${year}_${day}.nc HYCOM_GLBa0.08_${year}_${day}.nc
    rm -f HYCOM_GLBa0.08_ssh_${year}_${day}.nc HYCOM_GLBa0.08_temp_${year}_${day}.nc \ HYCOM_GLBa0.08_salt_${year}_${day}.nc HYCOM_GLBa0.08_u_${year}_${day}.nc \
HYCOM_GLBa0.08_v_${year}_${day}.nc
done
Gosh, should convert to bash.



Dear Kate

Thank you so much for sharing your works. I finished the 5th step and going to run make ic and bc file and got these error.

bc file
Code:
Traceback (most recent call last):
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/make_bdry_file.py", line 57, in <module>
    year = int(sys.argv[1])
IndexError: list index out of range
>>>


ic file (Ic problem has been solved).
Code:
Traceback (most recent call last):
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/make_ic_file.py", line 30, in <module>
    remap(file, 'temp', src_grd, dst_grd, dst_dir=dst_dir)
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/remap.py", line 127, in remap
    dxy=dxy, cdepth=cdepth, kk=kk)
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/Grid_HYCOM/flood_fast.py", line 107, in flood_fast
    varz[bottom[j,i]:,j,i] = varz[bottom[j,i],j,i]
IndexError: only integers, slices (`:`), ellipsis (`...`), numpy.newaxis (`None`) and integer or boolean arrays are valid indices
>>>


Thank you.
Regard.


Top
 Profile  
Reply with quote  
PostPosted: Sun Jun 23, 2019 4:53 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Quote:
year = int(sys.argv[1])
Did you give it a year as a command-line argument? It is expecting one. There's also a commented out line:
Code:
#lst_year = sys.argv[1:]
that can be enabled, allowing you to give it a list of years. Or you can edit the script to provide the year(s) internally.

For the other, I thought I had fixed all those "bottom" problems. Indeed, my copy of that line contains:
Code:
              varz[int(bottom[j,i]):,j,i] = varz[int(bottom[j,i]),j,i]


Top
 Profile  
Reply with quote  
PostPosted: Sun Jun 23, 2019 5:45 pm 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
Quote:
year = int(sys.argv[1])
Did you give it a year as a command-line argument? It is expecting one. There's also a commented out line:
Code:
#lst_year = sys.argv[1:]
that can be enabled, allowing you to give it a list of years. Or you can edit the script to provide the year(s) internally.

For the other, I thought I had fixed all those "bottom" problems. Indeed, my copy of that line contains:
Code:
              varz[int(bottom[j,i]):,j,i] = varz[int(bottom[j,i]),j,i]


Dear Kate
thank you for your answers, I put the year in then passed (when not python asked me for the year, I input but it still give me error)
the bc script run for a while then a message box appear to say that "Your program is still running! Do you want to kill it?" I said "No" for 4 time then the script seem to be stopped. May be I have some problems with these lines

Code:
processes = 4
p = Pool(processes)
# Trick to pass more than one arg
partial_do_file = partial(do_file, src_grd=src_grd, dst_grd=dst_grd)
results = p.map(partial_do_file, lst_file)


I have been run ROMS for multiple cores for several times. I have no idea about this.

Regard.


Top
 Profile  
Reply with quote  
PostPosted: Sun Jun 23, 2019 7:18 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Each system is different. I run on one with front-end nodes containing 8 cores, while the compute nodes contain 24-28 cores. I run these scripts on the front-end nodes because they have access to the filesystems containing these files, but I only ask for 4 cores in my "pool" so other people can get things done too. The limits on these two types of jobs are different, and you might have run into some system limits on cputime, etc. Do you have a sys admin to talk to about where and how to run jobs? These things can take time and my current ROMS-to-ROMS boundary extraction is taking about two days per year.


Top
 Profile  
Reply with quote  
PostPosted: Mon Jun 24, 2019 1:21 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
The limits on these two types of jobs are different, and you might have run into some system limits on cputime, etc. Do you have a sys admin to talk to about where and how to run jobs? .

Dear Kate
Thank you for help. I am trying to run all in my laptop before asking admin about my cputime. I run the script and have this message box



I notice that the cpu usage all are very low so may be the script was actually stopped.
please forgive me if I am a newbie in linux base systems so some time it makes me hard to understand the problems.
Have a nice week.
Regard


Attachments:
Screenshot from 2019-06-24 08-17-01.png
Screenshot from 2019-06-24 08-17-01.png [ 114.23 KiB | Viewed 670 times ]
Top
 Profile  
Reply with quote  
PostPosted: Mon Jun 24, 2019 5:18 am 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Maybe you should try a pool size of 1, not running in parallel. I'm not familiar with that exact display, but the swap seems high, like it doesn't all fit in memory and none of them are making any progress.


Top
 Profile  
Reply with quote  
PostPosted: Mon Jun 24, 2019 6:21 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
Maybe you should try a pool size of 1, not running in parallel. I'm not familiar with that exact display, but the swap seems high, like it doesn't all fit in memory and none of them are making any progress.


Dear Kate
Thanks for help. I put some flag into the code then

so I have error with p = Pool(processes), I also try with processes = 1 but the error stills. Could you please help me other way to skip using Pool?

Thank you.
Regard.


Attachments:
Screenshot from 2019-06-24 13-16-48.png
Screenshot from 2019-06-24 13-16-48.png [ 84.34 KiB | Viewed 631 times ]
Top
 Profile  
Reply with quote  
PostPosted: Mon Jun 24, 2019 6:34 am 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Does it object to processes = 1 the same way as processes = 4, or was there a different error? What if you don't let it kill the job?


Top
 Profile  
Reply with quote  
PostPosted: Mon Jun 24, 2019 7:46 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
Does it object to processes = 1 the same way as processes = 4, or was there a different error? What if you don't let it kill the job?


I have same err for 4 and 1 processes. If i dont kill it then the script hang like picture above with no err message and activity.


Top
 Profile  
Reply with quote  
PostPosted: Mon Jun 24, 2019 4:17 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Can you try on a system with more memory?


Top
 Profile  
Reply with quote  
PostPosted: Tue Jun 25, 2019 12:11 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
Can you try on a system with more memory?


I restart and try to run it and the err stills remain. My system is 8gb of ram with nearly 6gb free. May be problem is related to system open mpi libraries can not linked to anaconda packages? I have succeced running ROMS test cases with mpi run.
I try to run with this
Code:
print('start')
processes = 4
print('stop1')
p = Pool(processes)
# Trick to pass more than one arg
print('stop2')
partial_do_file = partial(do_file, src_grd=src_grd, dst_grd=dst_grd)
print('stop3')
results = p.map(partial_do_file, lst_file)
print('stop4')


when running with these stop1-4 flags, the message boxes appear 4 times, when stop1 printed I have 1st message box, I picked "No" 3 times then the stop2 printed out together with 4th message box, I continue pick "No" on 4th message box then the stop3 printed out then script hang and seem to be never printout stop4, at this time my system still have nearly 2gb ram free.

I changed the code from
Code:
partial_do_file = partial(do_file, src_grd=src_grd, dst_grd=dst_grd)
results = p.map(partial_do_file, lst_file)


to
Code:
results = p.map(partial(do_file(lst_file, src_grd=src_grd, dst_grd=dst_grd)))


then the script passed then I have this

Code:
Traceback (most recent call last):
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/make_bdry_file.py", line 101, in <module>
    results = p.map(partial(do_file(file=lst_file, src_grd=src_grd, dst_grd=dst_grd)))
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/make_bdry_file.py", line 30, in do_file
    zeta = remap_bdry(file, 'ssh', src_grd, dst_grd, dst_dir=dst_dir)
  File "/work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/remap_bdry.py", line 28, in remap_bdry
    dst_file = src_file.rsplit('/')[-1]
AttributeError: 'list' object has no attribute 'rsplit'


Regard.


Top
 Profile  
Reply with quote  
PostPosted: Wed Jun 26, 2019 12:36 am 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
You can read up on multiprocessing online. Something wasn't right about your last attempt.

I would check your system limits. With a bash shell:
Code:
(snowdrifts) [kshedstrom@chinook02 Tidal_bay]$ ulimit -a
core file size          (blocks, -c) 0
data seg size           (kbytes, -d) unlimited
scheduling priority             (-e) 0
file size               (blocks, -f) unlimited
pending signals                 (-i) 257144
max locked memory       (kbytes, -l) unlimited
max memory size         (kbytes, -m) unlimited
open files                      (-n) 1024
pipe size            (512 bytes, -p) 8
POSIX message queues     (bytes, -q) 819200
real-time priority              (-r) 0
stack size              (kbytes, -s) 10240
cpu time               (seconds, -t) unlimited
max user processes              (-u) 1024
virtual memory          (kbytes, -v) unlimited
file locks                      (-x) unlimited
and with the csh family:
Code:
chinook03.rcs.alaska.edu 318% limit
cputime      unlimited
filesize     unlimited
datasize     unlimited
stacksize    10240 kbytes
coredumpsize 0 kbytes
memoryuse    unlimited
vmemoryuse   unlimited
descriptors  1024
memorylocked unlimited
maxproc      1024


Honestly, what we really want is to have this using xarray and dask instead of multiprocessing. I don't (yet) know enough to do that rewrite.


Top
 Profile  
Reply with quote  
PostPosted: Wed Jun 26, 2019 1:47 am 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
You can read up on multiprocessing online. Something wasn't right about your last attempt.

I would check your system limits. With a bash shell:


My information

Code:
(pr2) bangdt@C2:~$ ulimit -a
core file size          (blocks, -c) 0
data seg size           (kbytes, -d) unlimited
scheduling priority             (-e) 0
file size               (blocks, -f) unlimited
pending signals                 (-i) 31039
max locked memory       (kbytes, -l) 16384
max memory size         (kbytes, -m) unlimited
open files                      (-n) 1024
pipe size            (512 bytes, -p) 8
POSIX message queues     (bytes, -q) 819200
real-time priority              (-r) 0
stack size              (kbytes, -s) 8192
cpu time               (seconds, -t) unlimited
max user processes              (-u) 31039
virtual memory          (kbytes, -v) unlimited
file locks                      (-x) unlimited

and
Code:
# limit
cputime    unlimited
filesize    unlimited
datasize    unlimited
stacksize    8192 kbytes
coredumpsize    0 kbytes
memoryuse    unlimited
memorylocked    16384 kbytes
maxproc    31039
openfiles    1024


Top
 Profile  
Reply with quote  
PostPosted: Wed Jun 26, 2019 4:39 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
Maybe this is biting you:
Code:
memorylocked    16384 kbytes


Top
 Profile  
Reply with quote  
PostPosted: Thu Jun 27, 2019 10:51 pm 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
Maybe this is biting you:
Code:
memorylocked    16384 kbytes


Dear Kate

by running this code

Code:
(pr2) root@C2:~# ulimit -l
16384
(pr2) root@C2:~# ulimit -l unlimited
(pr2) root@C2:~# ulimit -l
unlimited


then my system has
Code:
(pr2) root@C2:~# ulimit -a
core file size          (blocks, -c) 0
data seg size           (kbytes, -d) unlimited
scheduling priority             (-e) 0
file size               (blocks, -f) unlimited
pending signals                 (-i) 31040
max locked memory       (kbytes, -l) unlimited
max memory size         (kbytes, -m) unlimited
open files                      (-n) 1024
pipe size            (512 bytes, -p) 8
POSIX message queues     (bytes, -q) 819200
real-time priority              (-r) 0
stack size              (kbytes, -s) 8192
cpu time               (seconds, -t) unlimited
max user processes              (-u) 31040
virtual memory          (kbytes, -v) unlimited
file locks                      (-x) unlimited


and

Code:
(pr2) root@C2:~# csh
# limit
cputime    unlimited
filesize    unlimited
datasize    unlimited
stacksize    8192 kbytes
coredumpsize    0 kbytes
memoryuse    unlimited
memorylocked    unlimited
maxproc    31040
openfiles    1024


I assume that by unlock the memories size I could load the full lst_file in to the memory and run,
but when I try to run the script the problem stills remain.
may be stuck at this line and can not load the lst_file to call do_file.

Code:
results = p.map(partial_do_file, lst_file)


Regard


Top
 Profile  
Reply with quote  
PostPosted: Thu Jun 27, 2019 11:26 pm 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
What you want is to call do_file for each file in the list. Did it print the list of files to operate on? Is it what you expected? Then you can comment out the chunk with processes = 4 to the end. Add a loop going through the list of files, calling do_file(file, src_grd, dst_grd), something like:
Code:
for file in lst_file:
   do_file(file, src_grd, dst_grd)


Top
 Profile  
Reply with quote  
PostPosted: Sat Jun 29, 2019 1:17 pm 
Offline

Joined: Sun Mar 10, 2019 2:34 pm
Posts: 28
Location: University of Tsukuba
kate wrote:
What you want is to call do_file for each file in the list. Did it print the list of files to operate on? Is it what you expected? Then you can comment out the chunk with processes = 4 to the end. Add a loop going through the list of files, calling do_file(file, src_grd, dst_grd), something like:
Code:
for file in lst_file:
   do_file(file, src_grd, dst_grd)


Dear Kate
Thank you so much for your suggestion. Now I have all ic and bc created. But I still have problem with them.

1- All my ic and bc file are 2d grid file, not 3d comparing to the source data
2- My bc file doesnt have temp, salt, u and v data, only have 2d grid data in mask on v and u, and all the mask on u, v are the same value (=1), the bathymetry data in the bc file are ok.

my grid id is

Code:
# gridid definition file
#
id      = EASTSEA
name    = EASTSEA
grdfile = /work/apps/anaconda3/envs/pr2/lib/python2.7/site-packages/pyroms_toolbox/works/East-Sea_HYCOM/East_sea_400.nc
N       = 20
grdtype = roms
Vtrans  = 1
theta_s = 5
theta_b = 0.4
Tcline  = 3


I am not sure if my ic and bc file are ok or how to check it (I use IDV to view them)
https://drive.google.com/drive/folders/1Uk3OeHw1CtWaVGjJLZx5nz_quKu7JUq-?usp=sharing
Regard.


Top
 Profile  
Reply with quote  
PostPosted: Sun Jun 30, 2019 1:12 am 
Offline
User avatar

Joined: Wed Jul 02, 2003 5:29 pm
Posts: 3633
Location: IMS/UAF, USA
I like to use ncview to get a quick and dirty view of what's in a netCDF file. Looking at 'ncdump -h' on a file can also be useful. If you didn't get temp and salt, then something didn't go right in the boundary extraction.

If you like the idea of using these Python scripts, perhaps it's time to learn what you can about Python/numpy and the rest. Check out pdb, the Python debugger.


Top
 Profile  
Reply with quote  
Display posts from previous:  Sort by  
Post new topic Reply to topic  [ 32 posts ] 

All times are UTC


Who is online

Users browsing this forum: No registered users and 2 guests


You cannot post new topics in this forum
You cannot reply to topics in this forum
You cannot edit your posts in this forum
You cannot delete your posts in this forum
You cannot post attachments in this forum

Search for:
Jump to:  
Powered by phpBB® Forum Software © phpBB Group