NOAA-ORR-ERD / PyGnome

The General NOAA Operational Modeling Environment
https://gnome.orr.noaa.gov/doc/pygnome/index.html
Other
58 stars 44 forks source link

Gnome output with time sequence #27

Closed dkobashi closed 7 years ago

dkobashi commented 7 years ago

Hi,

I would like to run PyGnome using ROMS output and when I tried to run PyGnome for non-weathering particles with 5 day forecast model outputs, the dimension of the locations (longitude and latitude) did not have time dimension. The metadata of the gnome netcdf output is as follow. I was expecting to get the dimension of the particle locations as (time,data) rather than (data). And the same for depth. What did I do wrong? I checked script_ny_roms script to create my own run script. It created multiple figure files, though.

Thanks,

DJ

DJ Kobashi Dept. of Oceanography Texas A&M University

netcdf txla_gnome { dimensions: time = UNLIMITED ; // (121 currently) data = UNLIMITED ; // (3441 currently) two = 2 ; three = 3 ; variables: double time(time) ; time:comment = "unspecified time zone" ; time:long_name = "time since the beginning of the simulation" ; time:standard_name = "time" ; time:calendar = "gregorian" ; time:units = "seconds since 2017-04-10T00:00:00" ; int particle_count(time) ; particle_count:units = "1" ; particle_count:long_name = "number of particles in a given timestep" ; particle_count:ragged_row_count = "particle count at nth timestep" ; ushort spill_num(data) ; spill_num:long_name = "spill to which the particle belongs" ; spill_num:spills_map = "0: Spill, " ; double longitude(data) ; longitude:units = "degrees_east" ; longitude:long_name = "longitude of the particle" ; longitude:standard_name = "longitude" ; int age(data) ; age:units = "seconds" ; age:long_name = "age of particle from time of release" ; double depth(data) ; depth:units = "meters" ; depth:long_name = "particle depth below sea surface" ; depth:standard_name = "depth" ; depth:axis = "z positive down" ; double mass(data) ; mass:units = "kilograms" ; mass:long_name = "mass of particle" ; double latitude(data) ; latitude:units = "degrees_north" ; latitude:long_name = "latitude of the particle" ; latitude:standard_name = "latitude" ; short status_codes(data) ; status_codes:flag_meanings = "0: not_released, 2: in_water, 3: on_land, 7: off_maps, 10: evaporated, 12: to_be_removed," ; status_codes:long_name = "particle status code" ; status_codes:flag_values = "7 12 0 10 2 3" ; uint id(data) ; id:long_name = "particle ID" ;

// global attributes: :comment = "Particle output from the NOAA PyGnome model" ; :feature_type = "particle_trajectory" ; :conventions = "CF-1.6" ; :source = "PyGnome version v0.0.4" ; :references = "TBD" ; :institution = "NOAA Emergency Response Division" ; :creation_date = "2017-04-11T09:50:27" ;

ChrisBarker-NOAA commented 7 years ago

GNOME outputs the "particle" data in the nc_particles format:

http://noaa-orr-erd.github.io/nc_particles/nc_particle_standard.html

We have a Python module for working with that format in teh GnomeTool respository, in the post_gnome pacakge:

https://github.com/NOAA-ORR-ERD/GnomeTools

Why all this complication?

In short -- GNOME, and other particle tracking models, create and remove particles during the run, and does not, in general, know how many particles might ultimately get created. So they need to be stored in a "ragged_array" format.

netcdf does not support such a format directly -- thus the complication.

dkobashi commented 7 years ago

Dear Chris,

Thanks for your reply. My experience in particle tracking model is all about passive particles. So I did not know the Gnome model "creates/removes" particles during the run. Do you have some references that explain the mechanism? I would like to know more about it. As far as I checked, the Gnome manual does not articulate it (or is it?)

Ok. I've already installed GnomeTools and checked HowTo for nc_particles by checking the Boston example. When I tried to process my Gnome netcdf output, it did not work. When I tried to read position and trajectory data for a certain time step as with the case for the example, I just got an empty array.

import nc_particles idx=20 r = nc_particles.Reader('txla_gnome.nc') positions = r.get_timestep(idx) trajectory = r.get_individual_trajectory(idx)

The above line gave an empty array. The ROMS output is 1 hourly for 5 simulation days, so there should be data at time index 20.

I attach my run script below. I wonder if you could point out what the problem is for my PyGNome run. What I did was to load coastline file and load ROMS netcdf file via OPeNDAP. Then I ran PyGnome. No wind input. Many thanks in advance. -DJ

import os from datetime import datetime, timedelta

import numpy as np

from gnome import scripting from gnome import utilities from gnome.basic_types import datetime_value_2d, numerical_methods

from gnome.utilities.remote_data import get_datafile

from gnome.model import Model

from gnome.map import MapFromBNA from gnome.environment import GridWind from gnome.environment import GridCurrent from gnome.spill import point_line_release_spill, SpatialRelease, Spill from gnome.movers import RandomMover, constant_wind_mover, GridCurrentMover

from gnome.outputters import Renderer from gnome.outputters import NetCDFOutput from gnome.movers.py_current_movers import PyGridCurrentMover from gnome.movers.py_wind_movers import PyWindMover

import gnome.utilities.profiledeco as pd base_dir = os.path.dirname(file)

def make_model(images_dir=os.path.join(base_dir, 'images')): print 'initializing the model'

start_time = datetime(2017, 4, 10, 0, 0)

model = Model(start_time=start_time,
              duration=timedelta(days=5),
              time_step=timedelta(minutes=10))

mapfile = get_datafile(os.path.join(base_dir, 'txla_shelf.bna'))

print 'adding the map'
'''TODO: sort out MapFromBna's map_bounds parameter...
it does nothing right now, and the spill is out of bounds'''
model.map = MapFromBNA(mapfile, refloat_halflife=0.0)  # seconds

netcdf_file = 'txla_gnome.nc'
scripting.remove_netcdf(netcdf_file)
nc_outputter  = NetCDFOutput(netcdf_file, which_data='most',output_timestep=timedelta(hours=1))
print 'adding outputters'
model.outputters += nc_outputter

print 'adding a spill'
# for now subsurface spill stays on initial layer
# - will need diffusion and rise velocity
# - wind doesn't act
# - start_position = (-76.126872, 37.680952, 5.0),
particle_type='passive'
if particle_type=='passive':
    import ltrans_sub
    grdfile='/home/dj/Modeling/ROMS/ROMStools/mypython/grd/txla2_grd_v4_test_lcut_hglo_wtype.nc'
    get_locations=False
    save_file=True
    inpdir='./'
    outfile='txla_locs.txt'
    if get_locations:
        xx,yy,zz=ltrans_sub.create_locations_box(x=None,y=None,z=0.0,b=0.0,
                                          grdfile=grdfile,inpdir=inpdir,outfile=outfile,dx=10000.0,save_file=save_file)
    else:
        data=np.loadtxt(inpdir+outfile,delimiter=',')
        xx=data[:,0];yy=data[:,1];zz=data[:,2]
    positions=(xx,yy,zz)
    release = SpatialRelease(release_time=start_time,
                             start_position=positions)
    spill1  = Spill(release)
elif particle_type=='oil':
    spill1 = point_line_release_spill(num_elements=1000,
                                  start_position=(-94.5,29.0,0.0),
                                  release_time=start_time)

model.spills += spill1

print 'adding a RandomMover:'
model.movers += RandomMover(diffusion_coef=10000)

print 'adding a current mover:'
url = ('http://barataria.tamu.edu:8080/thredds/dodsC/NcML/oof_latest_agg')
cf = GridCurrent.from_netCDF(url)
u_mover = PyGridCurrentMover(cf, default_num_method='Euler')
model.movers += u_mover

return model

if name == "main": startTime = datetime.now() model = make_model() print "doing full run" model.full_run() print datetime.now() - startTime

ChrisBarker-NOAA commented 7 years ago

So I did not know the Gnome model "creates/removes" particles during the

run. Do you have some references that explain the mechanism? I would like to know more about it. As far as I checked, the Gnome manual does not articulate it (or is it?)

It may not be clear in there, but currently, GNOME adds particles as they are released, and removes then when they go off the map. So if you have an instantaneous release, and no particles go off the map, then the particle number will be constant.

Ok. I've already installed GnomeTools and checked HowTo for nc_particles by checking the Boston example. When I tried to process my Gnome netcdf output, it did not work. When I tried to read position and trajectory data for a certain time step as with the case for the example, I just got an empty array.

import nc_particles idx=20 r = nc_particles.Reader('txla_gnome.nc') positions = r.get_timestep(idx)

This should give you the positions of the particles at timestep 20

what do you get?

trajectory = r.get_individual_trajectory(idx)

this gives you the positions of particle number 20 -- do you ave more than 20 particles in your simulation?

The above line gave an empty array. The ROMS output is 1 hourly for 5 simulation days, so there should be data at time index 20.

you can try:

r.times

to be sure.

If the spill is out of bounds, you may not get any output!

-CHB

I attach my run script below. I wonder if you could point out what the

problem is for my PyGNome run. What I did was to load coastline file and load ROMS netcdf file via OPeNDAP. Then I ran PyGnome. No wind input. Many thanks in advance. -DJ

import os from datetime import datetime, timedelta

import numpy as np

from gnome import scripting from gnome import utilities from gnome.basic_types import datetime_value_2d, numerical_methods

from gnome.utilities.remote_data import get_datafile

from gnome.model import Model

from gnome.map import MapFromBNA from gnome.environment import GridWind from gnome.environment import GridCurrent from gnome.spill import point_line_release_spill, SpatialRelease, Spill from gnome.movers import RandomMover, constant_wind_mover, GridCurrentMover

from gnome.outputters import Renderer from gnome.outputters import NetCDFOutput from gnome.movers.py_current_movers import PyGridCurrentMover from gnome.movers.py_wind_movers import PyWindMover

import gnome.utilities.profiledeco as pd base_dir = os.path.dirname(file)

def make_model(images_dir=os.path.join(base_dir, 'images')): print 'initializing the model'

start_time = datetime(2017, 4, 10, 0, 0)

model = Model(start_time=start_time, duration=timedelta(days=5), time_step=timedelta(minutes=10))

mapfile = get_datafile(os.path.join(base_dir, 'txla_shelf.bna'))

print 'adding the map' '''TODO: sort out MapFromBna's map_bounds parameter... it does nothing right now, and the spill is out of bounds''' model.map = MapFromBNA(mapfile, refloat_halflife=0.0) # seconds

netcdf_file = 'txla_gnome.nc' scripting.remove_netcdf(netcdf_file) nc_outputter = NetCDFOutput(netcdf_file, which_data='most',output_timestep=timedelta(hours=1)) print 'adding outputters' model.outputters += nc_outputter

print 'adding a spill'

for now subsurface spill stays on initial layer

- will need diffusion and rise velocity

- wind doesn't act

- start_position = (-76.126872, 37.680952, 5.0),

particle_type='passive' if particle_type=='passive': import ltrans_sub grdfile='/home/dj/Modeling/ROMS/ROMStools/mypython/grd/txla2_grd_v4_test_lcut_hglo_wtype.nc' get_locations=False save_file=True inpdir='./' outfile='txla_locs.txt' if get_locations: xx,yy,zz=ltrans_sub.create_locations_box(x=None,y=None,z=0.0,b=0.0, grdfile=grdfile,inpdir=inpdir,outfile=outfile,dx=10000.0,save_file=save_file) else: data=np.loadtxt(inpdir+outfile,delimiter=',') xx=data[:,0];yy=data[:,1];zz=data[:,2] positions=(xx,yy,zz) release = SpatialRelease(release_time=start_time, start_position=positions) spill1 = Spill(release) elif particle_type=='oil': spill1 = point_line_release_spill(num_elements=1000, start_position=(-94.5,29.0,0.0), release_time=start_time)

model.spills += spill1

print 'adding a RandomMover:' model.movers += RandomMover(diffusion_coef=10000)

print 'adding a current mover:' url = ('http://barataria.tamu.edu:8080/thredds/dodsC/NcML/oof_latest_agg') cf = GridCurrent.from_netCDF(url) u_mover = PyGridCurrentMover(cf, default_num_method='Euler') model.movers += u_mover

return model

if name == "main": startTime = datetime.now() model = make_model() print "doing full run" model.full_run() print datetime.now() - startTime

— You are receiving this because you commented. Reply to this email directly, view it on GitHub https://github.com/NOAA-ORR-ERD/PyGnome/issues/27#issuecomment-293346638, or mute the thread https://github.com/notifications/unsubscribe-auth/AA38YG3tqVZNo8D61v3JyLZkQL2C2cGkks5ru77xgaJpZM4M6Of5 .

--

Christopher Barker, Ph.D. Oceanographer

Emergency Response Division NOAA/NOS/OR&R (206) 526-6959 voice 7600 Sand Point Way NE (206) 526-6329 fax Seattle, WA 98115 (206) 526-6317 main reception

Chris.Barker@noaa.gov

dkobashi commented 7 years ago

Dear Chris, Thanks for your reply.

import nc_particles idx=20 r = nc_particles.Reader('txla_gnome.nc') positions = r.get_timestep(idx) This should give you the positions of the particles at timestep 20

what do you get?

print positions {u'latitude': array([], dtype=float64), u'longitude': array([], dtype=float64)}

r.times to be sure.

print r.times gives all time values in datetime format and print r.get_all_timesteps() gives a value at the first time stamp and the rest of arrays are empty (see the output at the bottom) The particles are put every 10 km over the model domain (the Texas-Louisiana Shelf). So I am sure that many of the particles are still in the domain. I did not change anything for ROMS output in NetCDF format. Am I supposed to modify something? A strange thing is that the file size of the GNOME netcdf file is just 56KB. Too small.

In [3]: print r.times [datetime.datetime(2017, 4, 10, 0, 0) datetime.datetime(2017, 4, 10, 1, 0) datetime.datetime(2017, 4, 10, 2, 0) datetime.datetime(2017, 4, 10, 3, 0) datetime.datetime(2017, 4, 10, 4, 0) datetime.datetime(2017, 4, 10, 5, 0) datetime.datetime(2017, 4, 10, 6, 0) datetime.datetime(2017, 4, 10, 7, 0) datetime.datetime(2017, 4, 10, 8, 0) datetime.datetime(2017, 4, 10, 9, 0) datetime.datetime(2017, 4, 10, 10, 0) datetime.datetime(2017, 4, 10, 11, 0) datetime.datetime(2017, 4, 10, 12, 0) datetime.datetime(2017, 4, 10, 13, 0) datetime.datetime(2017, 4, 10, 14, 0) datetime.datetime(2017, 4, 10, 15, 0) datetime.datetime(2017, 4, 10, 16, 0) datetime.datetime(2017, 4, 10, 17, 0) datetime.datetime(2017, 4, 10, 18, 0) datetime.datetime(2017, 4, 10, 19, 0) datetime.datetime(2017, 4, 10, 20, 0) datetime.datetime(2017, 4, 10, 21, 0) datetime.datetime(2017, 4, 10, 22, 0) datetime.datetime(2017, 4, 10, 23, 0) datetime.datetime(2017, 4, 11, 0, 0) datetime.datetime(2017, 4, 11, 1, 0) datetime.datetime(2017, 4, 11, 2, 0) datetime.datetime(2017, 4, 11, 3, 0) datetime.datetime(2017, 4, 11, 4, 0) datetime.datetime(2017, 4, 11, 5, 0) datetime.datetime(2017, 4, 11, 6, 0) datetime.datetime(2017, 4, 11, 7, 0) datetime.datetime(2017, 4, 11, 8, 0) datetime.datetime(2017, 4, 11, 9, 0) datetime.datetime(2017, 4, 11, 10, 0) datetime.datetime(2017, 4, 11, 11, 0) datetime.datetime(2017, 4, 11, 12, 0) datetime.datetime(2017, 4, 11, 13, 0) datetime.datetime(2017, 4, 11, 14, 0) datetime.datetime(2017, 4, 11, 15, 0) datetime.datetime(2017, 4, 11, 16, 0) datetime.datetime(2017, 4, 11, 17, 0) datetime.datetime(2017, 4, 11, 18, 0) datetime.datetime(2017, 4, 11, 19, 0) datetime.datetime(2017, 4, 11, 20, 0) datetime.datetime(2017, 4, 11, 21, 0) datetime.datetime(2017, 4, 11, 22, 0) datetime.datetime(2017, 4, 11, 23, 0) datetime.datetime(2017, 4, 12, 0, 0) datetime.datetime(2017, 4, 12, 1, 0) datetime.datetime(2017, 4, 12, 2, 0) datetime.datetime(2017, 4, 12, 3, 0) datetime.datetime(2017, 4, 12, 4, 0) datetime.datetime(2017, 4, 12, 5, 0) datetime.datetime(2017, 4, 12, 6, 0) datetime.datetime(2017, 4, 12, 7, 0) datetime.datetime(2017, 4, 12, 8, 0) datetime.datetime(2017, 4, 12, 9, 0) datetime.datetime(2017, 4, 12, 10, 0) datetime.datetime(2017, 4, 12, 11, 0) datetime.datetime(2017, 4, 12, 12, 0) datetime.datetime(2017, 4, 12, 13, 0) datetime.datetime(2017, 4, 12, 14, 0) datetime.datetime(2017, 4, 12, 15, 0) datetime.datetime(2017, 4, 12, 16, 0) datetime.datetime(2017, 4, 12, 17, 0) datetime.datetime(2017, 4, 12, 18, 0) datetime.datetime(2017, 4, 12, 19, 0) datetime.datetime(2017, 4, 12, 20, 0) datetime.datetime(2017, 4, 12, 21, 0) datetime.datetime(2017, 4, 12, 22, 0) datetime.datetime(2017, 4, 12, 23, 0) datetime.datetime(2017, 4, 13, 0, 0) datetime.datetime(2017, 4, 13, 1, 0) datetime.datetime(2017, 4, 13, 2, 0) datetime.datetime(2017, 4, 13, 3, 0) datetime.datetime(2017, 4, 13, 4, 0) datetime.datetime(2017, 4, 13, 5, 0) datetime.datetime(2017, 4, 13, 6, 0) datetime.datetime(2017, 4, 13, 7, 0) datetime.datetime(2017, 4, 13, 8, 0) datetime.datetime(2017, 4, 13, 9, 0) datetime.datetime(2017, 4, 13, 10, 0) datetime.datetime(2017, 4, 13, 11, 0) datetime.datetime(2017, 4, 13, 12, 0) datetime.datetime(2017, 4, 13, 13, 0) datetime.datetime(2017, 4, 13, 14, 0) datetime.datetime(2017, 4, 13, 15, 0) datetime.datetime(2017, 4, 13, 16, 0) datetime.datetime(2017, 4, 13, 17, 0) datetime.datetime(2017, 4, 13, 18, 0) datetime.datetime(2017, 4, 13, 19, 0) datetime.datetime(2017, 4, 13, 20, 0) datetime.datetime(2017, 4, 13, 21, 0) datetime.datetime(2017, 4, 13, 22, 0) datetime.datetime(2017, 4, 13, 23, 0) datetime.datetime(2017, 4, 14, 0, 0) datetime.datetime(2017, 4, 14, 1, 0) datetime.datetime(2017, 4, 14, 2, 0) datetime.datetime(2017, 4, 14, 3, 0) datetime.datetime(2017, 4, 14, 4, 0) datetime.datetime(2017, 4, 14, 5, 0) datetime.datetime(2017, 4, 14, 6, 0) datetime.datetime(2017, 4, 14, 7, 0) datetime.datetime(2017, 4, 14, 8, 0) datetime.datetime(2017, 4, 14, 9, 0) datetime.datetime(2017, 4, 14, 10, 0) datetime.datetime(2017, 4, 14, 11, 0) datetime.datetime(2017, 4, 14, 12, 0) datetime.datetime(2017, 4, 14, 13, 0) datetime.datetime(2017, 4, 14, 14, 0) datetime.datetime(2017, 4, 14, 15, 0) datetime.datetime(2017, 4, 14, 16, 0) datetime.datetime(2017, 4, 14, 17, 0) datetime.datetime(2017, 4, 14, 18, 0) datetime.datetime(2017, 4, 14, 19, 0) datetime.datetime(2017, 4, 14, 20, 0) datetime.datetime(2017, 4, 14, 21, 0) datetime.datetime(2017, 4, 14, 22, 0) datetime.datetime(2017, 4, 14, 23, 0) datetime.datetime(2017, 4, 15, 0, 0)]

and print r.get_all_timesteps()

In [8]: r.get_all_timesteps() Out[8]: {u'latitude': [array([-97.59047, -97.32068, -97.05088, ..., 0. , 0. , 0. ]), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64)], u'longitude': [array([-97.68041, -97.41061, -97.14081, ..., 0. , 0. , 0. ]), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64), array([], dtype=float64)]}

ChrisBarker-NOAA commented 7 years ago

you have a bunch of particles in the zeroth timestep, then none after that. I don't have time to look at your script right now, but maybe that's a hint?

ChrisBarker-NOAA commented 7 years ago

@jay-hennen: could you take a look at this guy's script and see if you see anything wrong?

@dkobashi: do you have a modest-sized current file you could send us?

note: is it possible that everything is going off the map right away? try running it with no map (that will get you a all-water, no land map that encompases the whole earth.

dkobashi commented 7 years ago

Yes. I only want to add particles once at the beginning of the model and then I want to track where the particles I put go with 3-D ocean currents as time goes by. Did I do something wrong?

dkobashi commented 7 years ago

@ChrisBarker-NOAA

You can get the data via OPeNDAP from the following. http://barataria.tamu.edu:8080/thredds/dodsC/NcML/oof_latest_agg

Or HTTP (direct download of netcdf file) The file size = 25 GB http://barataria.tamu.edu:8080/thredds/fileServer/NcML/oof_latest_agg

The file is a ROMS (Regional Ocean Modeling System) history file. I get the model output via OPeNDAP, I wonder if that could be the problem.

Thanks in advance.

DJ

dkobashi commented 7 years ago

I clipped the ROMS output file by extracting the first 5 time steps with NCO.

http://barataria.tamu.edu:8080/thredds/dodsC/oof_other/roms_his_f_latest_clipped.nc (OPeNDAP) http://barataria.tamu.edu:8080/thredds/fileServer/oof_other/roms_his_f_latest_clipped.nc (HTTP) total size is 920MB

ChrisBarker-NOAA commented 7 years ago

Yes. I only want to add particles once at the beginning of the model

I meant your output file has particles only at the zeroth time step. So the spill added them to the model, then removed them all right away.

-CHB

jay-hennen commented 7 years ago

I would double check the positions of your particles and make sure the arrays that they are what SpatialRelease expects. I noticed two lines in your output:

{u'latitude': [array([-97.59047, -97.32068, -97.05088, ..., 0. , 0. , 0. ]), u'longitude': [array([-97.68041, -97.41061, -97.14081, ..., 0. , 0. , 0. ]),

followed by empty arrays. To me this means that the particles were placed on the first step, but outside the map and therefore on the second step got removed. I'm not entirely sure what the output of r.get_all_timesteps() is quite supposed to look like, but i would be surprised if latitude and longitude were both -97.xx....

dkobashi commented 7 years ago

@jay-hennen. Good catch! I found the problem. The problem was start_position on SpatialRelease. I created x,y,z array as (3, 1000). PyGnome requires (1000,3) instead. I had to transpose the array.

Now PyGnome runs fine for surface particles, but not for the particles put on the bottom. When I ran PyGNOME with the particles put on the bottom. All particles quickly come to surface right after the model starts. I wonder why.

Initial particle locations

part_loc_bottom

Depth against time

depth_vs_time

Particle trajectory

part_snapshot_bottom

ChrisBarker-NOAA commented 7 years ago

On Tue, Apr 11, 2017 at 3:07 PM, jay-hennen notifications@github.com wrote:

{u'latitude': [array([-97.59047, -97.32068, -97.05088, ..., 0. , 0. , 0. ]),

u'longitude': [array([-97.68041, -97.41061, -97.14081, ..., 0. , 0. , 0. ]),

followed by empty arrays. To me this means that the particles were placed on the first step, but outside the map and therefore on the second step got removed. I'm not entirely sure what the output of r.get_all_timesteps() is quite supposed to look like, but i would be surprised if latitude and longitude were both -97.xx....

good eye Jay -- clearly a latitude of -97 is wrong :-)

-CHB

— You are receiving this because you were mentioned. Reply to this email directly, view it on GitHub https://github.com/NOAA-ORR-ERD/PyGnome/issues/27#issuecomment-293414605, or mute the thread https://github.com/notifications/unsubscribe-auth/AA38YNFdVyGrQwjamCRJFZlhFOkp1n-Iks5ru_mEgaJpZM4M6Of5 .

--

Christopher Barker, Ph.D. Oceanographer

Emergency Response Division NOAA/NOS/OR&R (206) 526-6959 voice 7600 Sand Point Way NE (206) 526-6329 fax Seattle, WA 98115 (206) 526-6317 main reception

Chris.Barker@noaa.gov