Page Menu
Home
c4science
Search
Configure Global Search
Log In
Files
F76065807
controller.py
No One
Temporary
Actions
Download File
Edit File
Delete File
View Transforms
Subscribe
Mute Notifications
Award Token
Subscribers
None
File Metadata
Details
File Info
Storage
Attached
Created
Tue, Aug 6, 00:58
Size
21 KB
Mime Type
text/x-python
Expires
Thu, Aug 8, 00:58 (1 d, 21 h)
Engine
blob
Format
Raw Data
Handle
19657845
Attached To
rCTRACKER ctracker3
controller.py
View Options
###########################################################################
# #
# Copyright 2017 Andrea Cimatoribus #
# EPFL ENAC IIE ECOL #
# GR A1 435 (Batiment GR) #
# Station 2 #
# CH-1015 Lausanne #
# Andrea.Cimatoribus@epfl.ch #
# #
# Alexandre Serex #
# alexandre.serex@epfl.ch #
# #
# This file is part of ctracker #
# #
# ctracker is free software: you can redistribute it and/or modify it #
# under the terms of the GNU General Public License as published by #
# the Free Software Foundation, either version 3 of the License, or #
# (at your option) any later version. #
# #
# ctracker is distributed in the hope that it will be useful, #
# but WITHOUT ANY WARRANTY; without even the implied warranty #
# of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. #
# See the GNU General Public License for more details. #
# #
# You should have received a copy of the GNU General Public License #
# along with ctracker. If not, see <http://www.gnu.org/licenses/>. #
# #
###########################################################################
import
configparser
import
importlib
import
json
from
time
import
sleep
from
collections
import
namedtuple
from
itertools
import
product
from
os.path
import
join
import
numpy
as
np
from
core.router
import
connect
from
mods.grid
import
Grid
from
mods.particleContainer
import
ParticleContainer
from
parameters.netCDFController
import
NetCDFController
,
dataset
from
parameters.threadedQueue
import
INWorker
,
OUTWorker
from
parameters.RandomGenerator
import
RandomGenerator
from
tools.grid_extractor
import
mitgcm
,
mitgcm_raw
class
Controller
:
def
__init__
(
self
,
config_file_name
):
self
.
config_file_name
=
config_file_name
self
.
config
=
{}
self
.
config_methods
=
{
"json"
:
self
.
config_json
,
"ini"
:
self
.
config_ini
}
# Load the configuration file
try
:
file_extension
=
self
.
config_file_name
.
split
(
'.'
)[
-
1
]
self
.
config_methods
.
get
(
file_extension
,
self
.
config_invalid
)()
except
:
print
(
"
\n
Failed to load configuration file."
)
raise
print
(
"
\n
Configuration loaded successfully."
)
# Initialize the grid
try
:
# Concatenates the common values and the grid only values in a dict
grid_data
=
dict
(
self
.
config
[
"Common"
],
**
self
.
config
[
"Grid"
])
grid_tuple
=
namedtuple
(
'GridTuple'
,
grid_data
.
keys
())(
**
grid_data
)
# extract a raw grid from mitgcm module
self
.
raw_grid
=
mitgcm_raw
(
grid_tuple
)
grid
=
mitgcm
(
grid_tuple
)
# concatenate mitgm data to common and grid configuration data
grid_data
=
dict
(
grid_data
,
**
grid
)
# Creates a custom namedtuple with the above dict and feeds it to the grid
self
.
grid
=
Grid
(
"grid"
,
namedtuple
(
'GridData'
,
grid_data
.
keys
())(
**
grid_data
))
except
:
print
(
"
\n
Failed at initializing grid."
)
raise
# Initialize the simulation data as a particle container
try
:
part_cont_data
=
dict
(
self
.
config
[
"Common"
],
**
self
.
config
[
"Sim"
])
self
.
_part_container
=
ParticleContainer
(
"part_container"
,
namedtuple
(
'ParticleContainerData'
,
part_cont_data
.
keys
())(
**
part_cont_data
),
self
.
grid
)
except
:
print
(
"
\n
Failed at initializing simulation data."
)
raise
mitgcm_data
=
dict
(
self
.
config
[
"Common"
],
**
self
.
config
[
"MITgcm"
])
self
.
_data
=
namedtuple
(
'ParticleContainerData'
,
mitgcm_data
.
keys
())(
**
mitgcm_data
)
# Initialize the input output necessities
try
:
self
.
init_mitgcm
()
except
:
print
(
"
\n
Failed at initializing controller."
)
raise
# Check that all parameters' configuration is valid
self
.
_parameters
=
[
self
.
grid
,
self
.
_part_container
]
try
:
self
.
is_config_valid
()
[
p
.
is_config_valid
()
for
p
in
self
.
_parameters
]
except
:
print
(
"
\n
Invalid configuration setup."
)
raise
print
(
"
\n
Initialization successful."
)
def
is_config_valid
(
self
):
if
self
.
gcm_geometry
not
in
(
"curvilinear"
,
"cartesian"
):
raise
ValueError
(
"Unrecognised MITgcm geometry"
)
if
(
self
.
out_dt
%
self
.
gcm_dt
)
!=
0
:
raise
ValueError
(
"The GCM output interval must be "
"a multiple of the model time step."
)
if
self
.
iters
.
size
<
2
:
raise
ValueError
(
"To start the computation, the simulation "
"must span in time at least two GCM outputs"
)
if
(
self
.
seed_interval
%
self
.
out_dt
)
!=
0
:
raise
ValueError
(
"The seeding interval must be "
"a multiple of the GCM step."
)
if
self
.
inds_seed
:
# not sure if that np.any(...) does anything really
if
np
.
any
(
self
.
x_seed
>
self
.
grid
.
imax
)
or
\
np
.
any
(
self
.
x_seed
<
0
)
or
\
np
.
any
(
self
.
y_seed
>
self
.
grid
.
jmax
)
or
\
np
.
any
(
self
.
y_seed
<
0
)
or
\
np
.
any
(
self
.
z_seed
>
self
.
grid
.
kmax
)
or
\
np
.
any
(
self
.
z_seed
<
0
):
raise
ValueError
(
"Seeding with indexes beyond grid limits"
)
if
self
.
outfreq
in
(
"always"
,
"cross"
):
raise
NotImplementedError
(
"Outputing data with always and cross not yet supported."
)
if
self
.
outfreq
not
in
(
"gcmstep"
,
"always"
,
"cross"
,
"custom"
):
raise
ValueError
(
"Outfreq not recognised"
)
if
self
.
infreq
not
in
(
"gcmstep"
,
"custom"
):
raise
ValueError
(
"Infreq not recognised"
)
return
True
def
init_mitgcm
(
self
):
'''
Initialize the gcm input / output
'''
self
.
gcm_directory
=
self
.
_data
.
gcm_directory
self
.
out_prec
=
self
.
_data
.
out_prec
self
.
gcm_geometry
=
self
.
_data
.
gcm_geometry
self
.
gcm_start
=
self
.
_data
.
gcm_start
self
.
gcm_dt
=
self
.
_data
.
gcm_dt
self
.
out_dt
=
self
.
_data
.
out_dt
self
.
subiters
=
int
(
self
.
_data
.
subiters
)
self
.
dstep
=
1.0
/
self
.
subiters
self
.
out_dt
=
self
.
_data
.
out_dt
self
.
dtmax
=
self
.
out_dt
*
self
.
dstep
self
.
u_root
=
join
(
self
.
gcm_directory
,
self
.
_data
.
u_root
)
self
.
v_root
=
join
(
self
.
gcm_directory
,
self
.
_data
.
v_root
)
self
.
w_root
=
join
(
self
.
gcm_directory
,
self
.
_data
.
w_root
)
self
.
is2D
=
self
.
_data
.
is2D
print
(
"
\n
Simulation in 2D : "
,
self
.
is2D
)
self
.
complevel
=
self
.
_data
.
complevel
print
(
"
\n
Identify seeding points"
)
self
.
inds_seed
=
self
.
_data
.
inds_seed
self
.
seed_interval
=
self
.
_data
.
seed_interval
print
(
"
\n
Seeding interval: {0}s"
.
format
(
self
.
seed_interval
))
# save to netcdf, this will be the output file of the simulation
# we require NETCDF4
outfile_name
=
self
.
_data
.
outfile
.
rstrip
(
".nc"
)
self
.
inoutfile
=
outfile_name
+
"_inout.nc"
self
.
runfile
=
outfile_name
+
"_run.nc"
xGp1
,
yGp1
=
self
.
_get_geometry
()
self
.
raw_grid
.
attrs
[
"MITgcm_dir"
]
=
self
.
gcm_directory
# we also store the full mesh with edges coordinates
self
.
raw_grid
.
coords
[
"i_p1"
]
=
(
"i_p1"
,
np
.
arange
(
xGp1
.
shape
[
1
]))
self
.
raw_grid
.
coords
[
"j_p1"
]
=
(
"j_p1"
,
np
.
arange
(
xGp1
.
shape
[
0
]))
self
.
raw_grid
.
coords
[
"XG_p1"
]
=
((
"j_p1"
,
"i_p1"
),
xGp1
)
self
.
raw_grid
.
coords
[
"YG_p1"
]
=
((
"j_p1"
,
"i_p1"
),
yGp1
)
self
.
raw_grid
.
to_netcdf
(
self
.
runfile
,
mode
=
"w"
,
format
=
"NETCDF4"
,
engine
=
"netcdf4"
)
self
.
gcm_start
=
np
.
datetime64
(
self
.
_data
.
gcm_start
)
print
(
"
\n
Reference date of GCM:"
,
self
.
gcm_start
)
self
.
start
=
np
.
datetime64
(
self
.
_data
.
start
)
print
(
"
\n
Begin particle tracking simulation at:"
,
self
.
start
)
self
.
end
=
np
.
datetime64
(
self
.
_data
.
end
)
print
(
"
\n
End particle tracking simulation at:"
,
self
.
end
)
self
.
seed_start
=
np
.
datetime64
(
self
.
_data
.
seed_start
)
print
(
"
\n
Start seeding particles at:"
,
self
.
seed_start
)
self
.
seed_end
=
np
.
datetime64
(
self
.
_data
.
seed_end
)
print
(
"
\n
Stop seeding particles at:"
,
self
.
seed_end
)
self
.
outfreq
=
self
.
_data
.
outfreq
self
.
out_gcmstep
=
self
.
_data
.
out_gcmstep
# determine size of output buffer
self
.
nvals
=
1
if
self
.
outfreq
==
"always"
:
self
.
nvals
=
10
*
self
.
subiters
elif
self
.
outfreq
==
"cross"
:
self
.
nvals
=
15
# this should be a safe maximum
# first GCM iteration
iter_i
=
np
.
floor
(
(
self
.
start
-
self
.
gcm_start
)
.
astype
(
"timedelta64[s]"
)
.
astype
(
float
)
/
self
.
gcm_dt
)
.
astype
(
"int32"
)
# last GCM iteration
iter_e
=
np
.
ceil
(
(
self
.
end
-
self
.
gcm_start
)
.
astype
(
"timedelta64[s]"
)
.
astype
(
float
)
/
self
.
gcm_dt
)
.
astype
(
"int32"
)
self
.
iters
=
np
.
arange
(
iter_i
,
iter_e
+
1
,
int
(
self
.
out_dt
/
self
.
gcm_dt
))
self
.
infreq
=
self
.
_data
.
infreq
# we seed along a transect off the Rhone outflow
name_in
=
self
.
_data
.
in_data
.
rstrip
(
".py"
)
# import seed points coordinates
in_data
=
getattr
(
importlib
.
__import__
(
name_in
,
fromlist
=
[
name_in
]),
name_in
)
self
.
ijk_indseed
,
self
.
ijk_seed
=
self
.
_ijkseed
(
in_data
.
x_seed
,
in_data
.
y_seed
,
in_data
.
z_seed
,
self
.
inds_seed
)
self
.
xyz_seed
=
np
.
zeros
(
self
.
ijk_seed
.
shape
)
self
.
_part_container
.
grid
.
transform
(
self
.
ijk_seed
[:,
0
],
self
.
ijk_seed
[:,
1
],
self
.
ijk_seed
[:,
2
],
self
.
xyz_seed
)
self
.
n_seed
=
self
.
ijk_seed
.
shape
[
0
]
print
(
"
\n
Number of seeding points: {0}"
.
format
(
self
.
n_seed
))
# seeding steps
seed_i
=
np
.
floor
(
(
self
.
seed_start
-
self
.
gcm_start
)
.
astype
(
"timedelta64[s]"
)
.
astype
(
float
)
/
self
.
gcm_dt
)
.
astype
(
"int32"
)
seed_e
=
np
.
ceil
(
(
self
.
seed_end
-
self
.
gcm_start
)
.
astype
(
"timedelta64[s]"
)
.
astype
(
float
)
/
self
.
gcm_dt
)
.
astype
(
"int32"
)
iter_seed
=
np
.
arange
(
seed_i
,
seed_e
+
1
,
int
(
self
.
seed_interval
/
self
.
gcm_dt
))
self
.
timestamps_in
=
np
.
intersect1d
(
self
.
iters
,
iter_seed
)
print
(
"
\n
Number of seeding time steps:
%d
"
%
self
.
timestamps_in
.
size
)
self
.
ntrackmax
=
len
(
self
.
timestamps_in
)
*
self
.
n_seed
print
(
"
\n
Total number of particles to be released:
%d
"
%
self
.
ntrackmax
)
self
.
chunk_id
=
self
.
_data
.
chunk_id
self
.
chunk_time
=
self
.
_data
.
chunk_time
self
.
timestamps_out
=
[]
if
self
.
outfreq
==
"gcmstep"
:
pass
elif
self
.
outfreq
==
"custom"
:
name_out
=
self
.
_data
.
out_custom
.
rstrip
(
".py"
)
self
.
timestamps_out
=
getattr
(
importlib
.
__import__
(
name_out
,
fromlist
=
[
name_out
]),
name_out
)
def
_ijkseed
(
self
,
ii
,
jj
,
kk
,
inds
):
'''
Define the ii, jj, kk indices from which particles will be
released.
ii, jj, kk: list, array, or scalar with indices (i, j, k) from which
particles will be released (must be integers).
inds: if True (default), ii, jj and kk are interpreted
as indices, otherwise they are interpreted as lists of
exact release positions'''
ii
=
np
.
atleast_1d
(
np
.
squeeze
(
ii
))
jj
=
np
.
atleast_1d
(
np
.
squeeze
(
jj
))
kk
=
np
.
atleast_1d
(
np
.
squeeze
(
kk
))
if
(
ii
.
size
!=
jj
.
size
)
or
(
ii
.
size
!=
kk
.
size
):
raise
ValueError
(
"ii, jj and kk must have all the same dimension."
)
if
inds
:
return
np
.
atleast_2d
(
np
.
asarray
(
zip
(
np
.
int16
(
ii
),
np
.
int16
(
jj
),
np
.
int16
(
kk
)))),
\
np
.
atleast_2d
(
np
.
asarray
(
zip
(
ii
,
jj
,
kk
)))
else
:
# if MITgcm coordinates are passed, we have to load the model grid
# and then translate into the "normalised" index coordinates of
# tracmass
def
xy2grid
(
x
,
y
,
dX
,
dY
,
C
,
S
):
nx
=
(
C
*
x
+
S
*
y
)
/
dX
ny
=
(
-
S
*
x
+
C
*
y
)
/
dY
return
nx
,
ny
xG
=
self
.
_part_container
.
grid
.
xG
yG
=
self
.
_part_container
.
grid
.
yG
dX
=
self
.
_part_container
.
grid
.
dX
dY
=
self
.
_part_container
.
grid
.
dY
if
self
.
gcm_geometry
==
"curvilinear"
:
cs
=
self
.
_part_container
.
grid
.
CS
sn
=
self
.
_part_container
.
grid
.
SN
zG
=
self
.
_part_container
.
grid
.
Z
iout
=
np
.
zeros
(
ii
.
size
)
*
np
.
nan
jout
=
np
.
zeros
(
ii
.
size
)
*
np
.
nan
kout
=
np
.
zeros
(
ii
.
size
)
*
np
.
nan
iindout
=
np
.
zeros
(
ii
.
size
,
dtype
=
"int16"
)
jindout
=
np
.
zeros
(
ii
.
size
,
dtype
=
"int16"
)
kindout
=
np
.
zeros
(
ii
.
size
,
dtype
=
"int16"
)
for
nn
,
(
xx
,
yy
,
zz
)
in
enumerate
(
zip
(
ii
,
jj
,
kk
)):
jseed
,
iseed
=
np
.
unravel_index
(
np
.
argmin
((
xx
-
xG
)
**
2
+
(
yy
-
yG
)
**
2
),
xG
.
shape
)
if
self
.
gcm_geometry
==
"curvilinear"
:
nx
,
ny
=
xy2grid
(
xx
-
xG
[
jseed
,
iseed
],
yy
-
yG
[
jseed
,
iseed
],
dX
[
jseed
,
iseed
],
dY
[
jseed
,
iseed
],
cs
[
jseed
,
iseed
],
sn
[
jseed
,
iseed
])
else
:
nx
=
(
xx
-
xG
[
jseed
,
iseed
])
/
dX
[
jseed
,
iseed
]
ny
=
(
yy
-
yG
[
jseed
,
iseed
])
/
dY
[
jseed
,
iseed
]
# check the computed nx and ny, and in case seek in
# nearby cell
if
(
nx
<
0
)
or
(
nx
>
1
)
or
(
ny
<
0
)
or
(
ny
>
1
):
for
ijtry
in
product
(
range
(
-
5
,
6
),
range
(
-
5
,
6
)):
inew
=
iseed
+
ijtry
[
0
]
jnew
=
jseed
+
ijtry
[
1
]
if
(
inew
<
0
)
or
(
jnew
<
0
)
\
or
(
inew
>=
self
.
_part_container
.
grid
.
imax
)
or
(
jnew
>=
self
.
_part_container
.
grid
.
jmax
):
continue
if
self
.
gcm_geometry
==
"curvilinear"
:
nx
,
ny
=
xy2grid
(
xx
-
xG
[
jnew
,
inew
],
yy
-
yG
[
jnew
,
inew
],
dX
[
jnew
,
inew
],
dY
[
jnew
,
inew
],
cs
[
jnew
,
inew
],
sn
[
jnew
,
inew
])
else
:
nx
=
(
xx
-
xG
[
jnew
,
inew
])
/
dX
[
jnew
,
inew
]
ny
=
(
yy
-
yG
[
jnew
,
inew
])
/
dY
[
jnew
,
inew
]
if
(
nx
>=
0
)
and
(
nx
<
1
)
and
(
ny
>=
0
)
and
(
ny
<
1
):
iseed
=
inew
jseed
=
jnew
break
else
:
raise
ValueError
(
"Could not find the point "
"(x=
%f
, y=
%f
)"
%
(
xx
,
yy
))
z_here
=
zG
[:,
jseed
,
iseed
]
if
(
zz
>
z_here
.
max
())
or
(
zz
<=
z_here
.
min
()):
print
(
"
\n
Point outside vertical bounds at "
"x,y,z=
%.2f
,
%.2f
,
%.2f
"
%
(
xx
,
yy
,
zz
))
continue
kk
=
np
.
where
(
z_here
>
zz
)[
0
][
-
1
]
nz
=
(
zz
-
z_here
[
kk
])
/
(
z_here
[
kk
+
1
]
-
z_here
[
kk
])
iindout
[
nn
]
=
iseed
jindout
[
nn
]
=
jseed
kindout
[
nn
]
=
kk
iout
[
nn
]
=
iseed
+
nx
jout
[
nn
]
=
jseed
+
ny
kout
[
nn
]
=
kk
+
nz
good
=
np
.
isfinite
(
iout
)
iindout
=
iindout
[
good
]
jindout
=
jindout
[
good
]
kindout
=
kindout
[
good
]
iout
=
iout
[
good
]
jout
=
jout
[
good
]
kout
=
kout
[
good
]
return
np
.
atleast_2d
(
np
.
asarray
(
list
(
zip
(
iindout
[
np
.
isfinite
(
iout
)],
jindout
[
np
.
isfinite
(
jout
)],
kindout
[
np
.
isfinite
(
kout
)])))),
\
np
.
atleast_2d
(
np
.
asarray
(
list
(
zip
(
iout
[
np
.
isfinite
(
iout
)],
jout
[
np
.
isfinite
(
jout
)],
kout
[
np
.
isfinite
(
kout
)]))))
# end _ijkseed
def
config_invalid
(
self
,
*
args
,
**
kwargs
):
raise
EnvironmentError
(
"Unsupported configuration file extension."
)
def
config_json
(
self
,
*
args
,
**
kwargs
):
with
open
(
self
.
config_file_name
)
as
f
:
self
.
config
=
json
.
load
(
f
)
def
config_ini
(
self
,
*
args
,
**
kwargs
):
config
=
configparser
.
ConfigParser
()
config
.
read
(
self
.
config_file_name
)
# Store in self.config the values of the config file.
# Example:
# [DEFAULT] hello = world
# will translate to :
# >>> self.config["DEFAULT"]["hello"]
# world
for
section_name
in
config
.
sections
():
section
=
config
[
section_name
]
self
.
config
[
section_name
]
=
{}
for
key
in
section
:
self
.
config
[
section_name
][
key
]
=
config
[
section_name
][
key
]
def
run
(
self
):
"""
Main function
"""
print
(
"
\n
Starting ctracker"
)
# open file for output
with
dataset
(
self
.
inoutfile
,
mode
=
"w"
)
as
nc_inout
,
\
dataset
(
self
.
runfile
,
mode
=
"r+"
)
as
nc_run
:
to_netcdf
=
NetCDFController
(
"tonetcdf"
,
nc_run
,
nc_inout
,
(
self
.
_part_container
.
iendw
,
self
.
_part_container
.
iende
,
self
.
_part_container
.
jendn
,
self
.
_part_container
.
jends
),
self
.
ntrackmax
,
self
.
gcm_start
,
self
.
complevel
,
self
.
chunk_id
,
self
.
chunk_time
)
in_worker
=
INWorker
(
"inworker"
,
self
.
iters
,
self
.
grid
,
(
self
.
u_root
,
self
.
v_root
,
self
.
w_root
),
self
.
is2D
,
self
.
out_prec
,
self
.
_part_container
.
ff
)
out_worker
=
OUTWorker
(
"outworker"
,
self
.
ijk_seed
,
self
.
xyz_seed
,
self
.
outfreq
)
print
(
"
\n
Start main loop
\n
"
)
seed_data
=
{
'nvals'
:
self
.
nvals
,
'n_seed'
:
self
.
n_seed
,
'ijk_seed'
:
self
.
ijk_seed
,
'ijk_indseed'
:
self
.
ijk_indseed
}
seed_ctx
=
namedtuple
(
'SeedCtx'
,
seed_data
.
keys
())(
**
seed_data
)
connect
(
in_worker
,
self
.
_part_container
)
connect
(
self
.
_part_container
,
out_worker
)
connect
(
out_worker
,
to_netcdf
)
connect
(
out_worker
,
to_netcdf
,
name_out
=
"sync"
,
name_in
=
"sync"
)
in_worker
.
start
()
out_worker
.
start
()
self
.
_part_container
.
run
(
self
.
iters
,
self
.
timestamps_in
,
seed_ctx
)
"""
generator = RandomGenerator('randgen')
out_worker = OUTWorker("outworker", self.ijk_seed, self.xyz_seed, self.outfreq)
connect(generator, out_worker)
connect(out_worker, to_netcdf)
connect(out_worker, to_netcdf, "sync", "sync")
n_part = 200
for i in range(0, 30):
generator.to_out_worker(65665800.0 + i*1800.0, (n_part, 1, 4), (n_part, 1, 3), (n_part,), (n_part,))
out_worker.process_output()
"""
# end with statement
print
(
"
\n
Simulation ended."
)
# end run
def
_get_geometry
(
self
):
"""
Convenience function returning XG and YG of mitgcm
"""
xG
=
np
.
zeros
((
self
.
grid
.
jmax
+
1
,
self
.
grid
.
imax
+
1
))
yG
=
np
.
zeros
((
self
.
grid
.
jmax
+
1
,
self
.
grid
.
imax
+
1
))
xG
[:
-
1
,
:
-
1
]
=
self
.
grid
.
xG
yG
[:
-
1
,
:
-
1
]
=
self
.
grid
.
yG
dxG
=
self
.
grid
.
dX
dyG
=
self
.
grid
.
dY
if
self
.
gcm_geometry
==
"curvilinear"
:
cs
=
self
.
grid
.
CS
sn
=
self
.
grid
.
SN
# Fill (approximate) end points of the grid
xG
[:
-
1
,
-
1
]
=
xG
[:
-
1
,
-
2
]
+
dxG
[:,
-
1
]
*
cs
[:,
-
1
]
xG
[
-
1
,
:
-
1
]
=
xG
[
-
2
,
:
-
1
]
-
dyG
[
-
1
,
:]
*
sn
[
-
1
,
:]
# we lack the last metric at the NE corner, so we use the
# nearby metric
xG
[
-
1
,
-
1
]
=
xG
[
-
1
,
-
2
]
+
dxG
[
-
1
,
-
1
]
*
cs
[
-
1
,
-
1
]
yG
[
-
1
,
:
-
1
]
=
yG
[
-
2
,
:
-
1
]
+
dyG
[
-
1
,
:]
*
cs
[
-
1
,
:]
yG
[:
-
1
,
-
1
]
=
yG
[:
-
1
,
-
2
]
+
dxG
[:,
-
1
]
*
sn
[:,
-
1
]
yG
[
-
1
,
-
1
]
=
yG
[
-
2
,
-
1
]
+
dyG
[
-
1
,
-
1
]
*
cs
[
-
1
,
-
1
]
elif
self
.
gcm_geometry
==
"cartesian"
:
# Fill (approximate) end points of the grid
xG
[:
-
1
,
-
1
]
=
xG
[:
-
1
,
-
2
]
+
dxG
[:,
-
1
]
xG
[
-
1
,
:]
=
xG
[
-
2
,
:]
yG
[
-
1
,
:
-
1
]
=
yG
[
-
2
,
:
-
1
]
+
dyG
[
-
1
,
:]
yG
[:,
-
1
]
=
yG
[:,
-
2
]
else
:
raise
ValueError
(
"Grid geometry not recognised."
)
return
xG
,
yG
Event Timeline
Log In to Comment