Opened 13 years ago

Closed 13 years ago

Last modified 13 years ago

#4143 closed enhancement (fixed)

GCP Refinement

Reported by: goocreations Owned by: warmerdam
Priority: normal Milestone: 1.9.0
Component: default Version: svn-trunk
Severity: normal Keywords: gcp refine refinement


Here is a patch for gdalwarp to refine GCPs by using a -refine_gcps flag. GCPs will be detected as outliers by using a Least Squares Fit and not be used during warping.

Attachments (4)

gcp_refinement.diff (14.3 KB ) - added by goocreations 13 years ago.
SVN patch for the GCP refinement
gcp_refinement_22644.diff (24.7 KB ) - added by goocreations 13 years ago.
Updated GCP refinement
gcp_refinement_22652.diff (24.9 KB ) - added by goocreations 13 years ago.
Updated GCP refinement
gcp_refinement_22669.diff (2.4 KB ) - added by goocreations 13 years ago.

Download all attachments as: .zip

Change History (19)

by goocreations, 13 years ago

Attachment: gcp_refinement.diff added

SVN patch for the GCP refinement

comment:1 by Even Rouault, 13 years ago

I've reviewed the patch (but not tested). My general feeling is that it looks good. However here are my remarks :

  • in apps/gdalwarp.cpp, the help "[-refine_gcps]" should mention it needs 2 args. --> [-refine_gcps minimumGcps tolerance].
  • "The tolerance is passed to adjust when a GCP will be eliminated." --> it doesn't mention in which units tolerance should be expressed : pixels or target SRS units ?
  • I think that the -refince_gcps option only works with polynomial interpolation, not -tps option, right ? If so this should be documented.
  • the 2 new transformation options, REFINE_MINIMUM and REFINE_TOLERANCE, should be documented in GDALCreateGenImgProjTransformer2() (alg/gdaltransformer.cpp)
  • I'm wondering if the minimumGcps parameter should be validated. Particularly if someone passes a too small value and the interpolation doesn't work. There's perhaps validation at lower stages, but I think you should try with small values like 0,1,2 and see what happens.
  • GDALSerializeGCPTransformer() and GDALDeserializeGCPTransformer() should be updated with the 2 new params, so that gdalwarp -of VRT works. Make sure that the deserializer still works with old VRTs that wouldn't have the 2 new options.
  • it would be good to stick to the reversed polish notations used generally in GDAL for variable naming. E.g. minimumGcps --> nMinimumGcsp (n for integer) and tolerance -> dfTolerance (df for double float). At least in files where it is used in order to keep consistency. (you'll tell me that alg/gdal_crs.c is an exception to this convention... --> it is third-party code originally). Also, we tend to use 4 spaces for each indentation level.
  • in worst_outlier(), residuals array is allocated but never freed.
  • I think I got the gist of it, but a bit more comment in worst_outlier() and remove_outliers() (especially in the former) would be appreciated to ease the understanding the maths involved.

We'd appreciate if you could attach a sample dataset to the ticket and a gdalwarp command line example that illustrates the benefit of the patch. Ideally your dataset would be small and freely redistribuable so that we can add a new regression test to the GDAL autotest suite.

As far as the target version, this is really a new feature involving new API, options, parameters, etc.., so 1.8.1 isn't really appropriate. So make sure your patch is based against latest trunk (or at least that it applies without rejects on it). Thanks !

comment:2 by goocreations, 13 years ago

I've attached an example. I have 51 GCPs, one of these is misplaced. The data shows an extreme case. Using gdalwarp will result in a strangely warped image, but -refine_gcps will remove this misplaced GCP. I'm still busy with the comments/improvements as suggested above, but here is the test data. Usage:

gdalwarp -order 3 test.tif out.tif (will result in inaccurate warping) gdalwarp -order 3 -refine_gcps 32 1.9 test.tif out.tif (will remove the misplaced GCP and result in accurate warping)

comment:3 by goocreations, 13 years ago


comment:4 by goocreations, 13 years ago

Hi rouault, hope you get this message. regarding the minimumGcps that should be validated: Currently if the user specifies any number lower than needed by the interpolation, the process is started, but then the user is presented with an error: "Failed to compute GCP transform: Not enough points available" (as specified in gdal_crs.c). Is this enough, or should there be a more advanced technique to calculate this before the process is started? I'm currently changing the code so that 0 is the minimum and if minimumGcps is higher than the actual number of GCPs, minimumGcps is set to the number of GCPs.

comment:5 by Even Rouault, 13 years ago

I think the current error message is fine. I just wanted to make sure there's no crash in this situation.

by goocreations, 13 years ago

Attachment: gcp_refinement_22644.diff added

Updated GCP refinement

comment:6 by goocreations, 13 years ago

I've updated the patch, build and tested it against revision 22644. Please note that the parameters for -refine_gcps have changed: the tolerance (in SRS units) is now first and the minimum GCPs the second parameter. If the user omits minimum_gcps, the minimum_gcps according to the current model will be chosen. The minimum GCPs are limited to [0, GCPs available]. If too many GCPs were removed and the model cannot be constructed, a default GDAL error of too few points will be shown to the user. I've updated everything and improved the comments. I hope everything is in order. Please let me know if I forgot something.

comment:7 by Even Rouault, 13 years ago

Hum, are you sure that it works correctly if no explicit minimumGcps is passed ? In alg/gdaltransformer.cpp, I see at line 1108 that you create GDALCreateGCPRefineTransformer() only if nMinimumGcps >= 0. But if no explicit minimumGcps is passed, nMinimumGcps is set to -1 if REFINE_MINIMUM_GCPS is unset or set to default. I think you would need a bRefine flag in that method, that would be set to TRUE if REFINE_TOLERANCE isn't NULL for example.

I also think that the test if(dfCurrentDifference < 0) at line 1046 of alg/gdal_crs.c will always evaluate to false. dfCurrentDifference is set just at the line before to padfResiduals[nI] which is a square root, thus >= 0.

by goocreations, 13 years ago

Attachment: gcp_refinement_22652.diff added

Updated GCP refinement

comment:8 by goocreations, 13 years ago

Sorry about that. Hope everything is in order now.

comment:9 by Even Rouault, 13 years ago

Resolution: fixed
Status: newclosed

r22655 /trunk/ (6 files in 4 dirs): gdalwarp: add -refine_gcps option to discard outliers GCPs before warping (#4143)

I'v made the following changes w.r.t your latest patch :

  • small fix when -refine_gcps is on the last position of the command line and the minimumGcps is not provided
  • fix memory leak of removed psInfo->pasGCPList[nIndex].pszId and psInfo->pasGCPList[nIndex].pszInfo in remove_outliers()
  • simplify code of remove_outliers() to avoid the creation of a temporary array in the loop
  • I've also added a dump autotest to check this new feature

Due to my changes, please check that it still works however.

comment:10 by goocreations, 13 years ago

There is a major flaw in my code. I've updated my comments, but never my actual code. So here is a small patch, only changes 3 lines of code. Without this the refinement won't behave as expected.

Sorry about this.

by goocreations, 13 years ago

Attachment: gcp_refinement_22669.diff added

comment:11 by goocreations, 13 years ago

I've also changed the gdalwarp comment to indicate that minimum_gcps is optional

comment:12 by riaanvddool, 13 years ago

Resolution: fixed
Status: closedreopened

comment:13 by Even Rouault, 13 years ago

Resolution: fixed
Status: reopenedclosed

r22670 /trunk/gdal/ (alg/gdal_crs.c apps/gdalwarp.cpp): GCP refinement: compare directly against tolerance value (gcp_refinement_22669.diff, #4143)

comment:14 by riaanvddool, 13 years ago

Question asked on irc:

[10:34] <riaanvddool> hi guys, just wondering if there is any info on when 1.9.0 would be released? [10:34] <riaanvddool> i am specifically interested in the new gcp_refine functionality [10:35] <riaanvddool> i would like to write a magazine article for a local (ZA) magazine on our contribution to this functionality, but it doesn't make much sense until the functionality is available in a release [10:36] <riaanvddool> obviously the article would also provide good publicity for GDAL [10:36] <riaanvddool> also: if any other test cases etc are needed, would also like to know about it

comment:15 by Even Rouault, 13 years ago

GDAL releases occur generally on a yearly basis, in late December or January. If people want to try it before, they can for example try the Windows daily builds (the ones tagged -development) :

Note: See TracTickets for help on using tickets.