pism/pism — Issue #10: Make vel2tauc friendly to PISM workflow
Author: damaxwell
Created: 2011-11-04 22:06:42+00:00
State: open
Labels: task

Specification of `vel2tauc.py`:

(As revised following a November 4 2011 conversation with [Ed](https://github.com/bueler))

The job of `vel2tillphi.py` is to take a PISM model state, a set of 'observed ssa velocities', a scalar function `vel_misfit_weight` (the variable previously known as `range_l2_weight`), a previous guess for `tauc`, and produce an inferred value of `tauc`.

Variables for the model state are specified by either `-i in_file` or `-a in_file`; we'll call `in_file` the input file.

Output is stored in the output file. If `-a in_file` is specified, the input file is the output file. Otherwise, `-o out_file` may be used to indicate an output file. It's an error to use both `-a` and `-o`. If `-i` is used but not `-o`, then the output is placed in `vel2tillphi_input_file`

Variables for the 'observed ssa velocities' and misfit weight are specified in the inverse data file. It can be specified by `-inv_data foo.nc`, and is otherwise taken to be the input file.

The following variables are meaningful in the inverse data file:

* u/v_surface_observed -- surface velocities from measurement
* vel_surface_mask -- = 1 if velocity measurement is present
* u/v_ssa_observed -- explained below
* vel_misfit_weight -- a nodewise weight used in the L2 misfit functional
* tauc_prior -- explained below

If the variables `u/v_ssa_observed` are present, these are taken as input for inversion. Otherwise, the variables `u/v_surface_observed` must be present, and a solution of the SIA, `u/v_sia_observed` is constructed using data from the input file. We set `u/v_ssa_observed = u/v_surface_observed - u/v_sia_observed`. One can imagine in the future perhaps creating `u/v_sia_observed` from geometry different from the geometry in the input file; hence the appendage observed here.

I will be experimenting with using a `vel_misfit_weight` that is not specified on a per-node basis but on a per-element basis. For now, I'll use `vel_elementwise_misfit_weight` for this quantity.

The inversion run requires a starting point for `tauc`. The default is to take this variable from the input file. It may be useful to be able to restart an inversion from an aborted attempt. During each iteration of the inversion, the current value of tauc produced will be saved in the output file as `tauc_inv`. If `-inv_restart` is specified, then the starting `tauc` is taken to be `tauc_inv` in the input file. An explicit value of `tauc` can be specified in the inverse data file as `tauc_prior`, and use of this starting point is indicated via `-use_tauc_prior`.

The final value of `tauc` computed from inversion is saved in the output file as `tauc`. The starting point for `tauc` is saved in the output file as `tauc_prior`. We save the data we used from the inverse data file to the output file, including `u/v_sia_observed` if they were generated. If not running in append mode, we also write all PISM model state variables required for the SSA/SIA to the output file.

At the end of a successful inversion run that started from a tauc and producing "output_file", we should be able to replicate *exactly* the same inversion with the flags :-i/a output_file -use_tauc_prior. The output file will be managed in such a way that if there is a crash or an abort, the inversion can be continued from the output_file with with -i/a output_file -inv_restart.

I envision that standard use will be

vel2tauc.py -a input-file.nc -inv_data vel-obs.nc [other flags]

with `tauc` from the input data being used as the initial guess and consequently saved as tauc_prior post-inversion.

After a successful inversion, the output file will contain:

tauc: the tauc generated from inversion
tauc_prior: the tauc starting point for the inversion
tauc_inv: the most recent iteration of tauc produced during an inversion run