Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Temporal gridshifting #1015

Merged
merged 6 commits into from May 24, 2018
Merged

Temporal gridshifting #1015

merged 6 commits into from May 24, 2018

Conversation

@kbevers
Copy link
Member

@kbevers kbevers commented May 24, 2018

Extending the gridshift operations into the fourth dimension. The idea is to be able to use the grid shift operations as step functions when making transformations in areas subject to deformation after seismic events. This is useful when transforming coordinate from one epoch to another, e.g. in dynamic reference frames. This work is closely related to the proposed deformation model grid format in #1001.

I will expand the docs later to today, for now I am just checking that the code builds properly.

kbevers added 3 commits May 7, 2018

if (!pj_param(P->ctx, P->params, "tgrids").i) {
proj_log_error(P, "vgridshift: +grids parameter missing.");
return pj_default_destructor(P, PJD_ERR_NO_ARGS);
}

if (pj_param(P->ctx, P->params, "tt_final").i) {
Q->t_final = pj_param (P->ctx, P->params, "dt_final").f;

This comment has been minimized.

@rouault

rouault May 24, 2018
Member

we should likely have a pj_eval_param_time(const char*) function to do the below processing, since it is repeated above.
localtime() is known on POSIX to be thread-unsafe. It is better to use localtime_r() instead, when available

This comment has been minimized.

@kbevers

kbevers May 24, 2018
Author Member

we should likely have a pj_eval_param_time(const char*) function to do the below processing, since it is repeated above.

Good point. I intend to add similar functionality to deformation and possibly helmert as well. I will take care of it then.

localtime() is known on POSIX to be thread-unsafe. It is better to use localtime_r() instead, when available

When can localtime_r be expected to be available?

This comment has been minimized.

@rouault

rouault May 24, 2018
Member

When can localtime_r be expected to be available?

We'd need a configure / cmake check to test for its availability. Or alternatively "port" GDAL's https://github.com/OSGeo/gdal/blob/master/gdal/port/cpl_time.cpp#L64 which as indicated in the file header comes from public domain code.

This comment has been minimized.

@kbevers

kbevers May 24, 2018
Author Member

I think it is a bit out of scope for this PR. How big a risk is it to include the code as it is now? localtime() is only called once during creation of the PJ object.

This comment has been minimized.

@rouault

rouault May 24, 2018
Member

The risk is probably low. A TODO in the code + issue to remember should be good enough

This comment has been minimized.

@kbevers

kbevers May 24, 2018
Author Member

Created #1016

This comment has been minimized.

@busstoptaktik

busstoptaktik May 24, 2018
Member

You can probably significantly reduce (but not eliminate) the risk by saying:

struct tm cucumber;
...
cucumber = *localtime(...);

Essentially doing a copy of the the localtimeinternal struct tm, rather than storing a pointer to it.

struct pj_opaque_hgridshift *Q = (struct pj_opaque_hgridshift *) P->opaque;
PJ_COORD point = obs;
if (Q->t_final != 0.0 && Q->t_epoch != 0.0) {
if (obs.lpzt.t < Q->t_epoch && Q->t_final > Q->t_epoch)

This comment has been minimized.

@rouault

rouault May 24, 2018
Member

shouldn't that be be <= and >= ?

Or so if t_final and t_epoch are specified, we apply the grid shift only if the observation time is below the epoch, and the epoch below the final epoch. That doesn't immediately seem obvious to me. Should probably be explained in the .rst files

This comment has been minimized.

@kbevers

kbevers May 24, 2018
Author Member

I have just updated the .rst files. Do they cover the topic better now?

From a physical point of view I would say that the case where obs.lpzt.t == Q->t_epoch is undefined, since the shift doesn't happen instantaneously but rather over a timespan of minutes to days. To me it makes the most sense to say that the offset is applied only if the observation epoch is before the event that caused the deformation. Should you happen to have acquired a coordinate during the deformation event it will be partly to fully shifted and we dont' really have the means to do anything about that.

A similar point can be made for the second part of the condition, but I think from a practical point of view it make sense to change it to Q->t_final >= Q->t_epoch.

This comment has been minimized.

@rouault

rouault May 24, 2018
Member

Yes the docs make things really clear (sorry I started looking at the code before reading your PR general comment where you mentioned that would be later detailed :-) )
And actually no strong opinion if that should be strict or equal comparisons.

@@ -11,9 +11,11 @@ Change Vertical datum change by grid shift
+-----------------+--------------------------------------------------------------------+
| **Domain** | 3D |
+-----------------+--------------------------------------------------------------------+
| **Input type** | Geodetic coordinates (horizontal), meters (vertical). |
| **Input type** | Geodetic coordinates (horizontal), meters (vertical), |

This comment has been minimized.

@rouault

rouault May 24, 2018
Member

hgridshift.rst should also be similarly updated

Copy link
Member

@busstoptaktik busstoptaktik left a comment

This review consists of just one comment about code flow and readability, but it applies 4 places :-)

return obs;
static PJ_COORD reverse_4d(PJ_COORD obs, PJ *P) {
struct pj_opaque_hgridshift *Q = (struct pj_opaque_hgridshift *) P->opaque;
PJ_COORD point = obs;

This comment has been minimized.

@busstoptaktik

busstoptaktik May 24, 2018
Member

I think the code below (and the similar in the forward_4d case) will be more readable if the conditionals are inverted, so the if/if/else-branches can be decoupled:

/* If transformation is not time restricted, we always call it */
if (Q->t_final==0 || Q->t_epoch==0) {
    point.lpz = reverse_3d (obs.xyz, P);
    return point;
}

/* Time restricted - only apply transform if within time bracket */
if (obs.lpzt.t < Q->t_epoch && Q->t_final > Q->t_epoch)
    point.lpz = reverse_3d (obs.xyz, P);
return point;

This comment has been minimized.

@kbevers

kbevers May 24, 2018
Author Member

Good suggestion. Updated the code.

PJ_COORD point = {{0,0,0,0}};
point.xyz = forward_3d (obs.lpz, P);
struct pj_opaque_vgridshift *Q = (struct pj_opaque_vgridshift *) P->opaque;
PJ_COORD point = obs;

This comment has been minimized.

@busstoptaktik

busstoptaktik May 24, 2018
Member

cf. comment in PJ_hgridshift.c regarding inverted conditionals

PJ_COORD point;
point.lpz = reverse_3d (obs.xyz, P);
struct pj_opaque_vgridshift *Q = (struct pj_opaque_vgridshift *) P->opaque;
PJ_COORD point = obs;

This comment has been minimized.

@busstoptaktik

busstoptaktik May 24, 2018
Member

cf. comment in PJ_hgridshift.c regarding inverted conditionals

@kbevers kbevers mentioned this pull request May 24, 2018
0 of 3 tasks complete
@kbevers kbevers force-pushed the kbevers:temporal-gridshift branch from 9322ec2 to 1e03829 May 24, 2018
@kbevers
Copy link
Member Author

@kbevers kbevers commented May 24, 2018

@rouault and @busstoptaktik, thanks for the swift and helpful reviews. I am merging this now so that it can be included in 5.1.0.

@kbevers kbevers merged commit 0b5c362 into OSGeo:master May 24, 2018
3 checks passed
3 checks passed
continuous-integration/appveyor/pr AppVeyor build succeeded
Details
continuous-integration/travis-ci/pr The Travis CI build passed
Details
coverage/coveralls Coverage increased (+0.06%) to 77.649%
Details
@kbevers kbevers added this to the 5.1.0 milestone May 24, 2018
@kbevers kbevers deleted the kbevers:temporal-gridshift branch Jul 9, 2018
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Projects
None yet
Linked issues

Successfully merging this pull request may close these issues.

None yet

3 participants
You can’t perform that action at this time.