[GRASS-SVN] r44653 - in
grass/branches/releasebranch_6_4/imagery/i.ortho.photo: .
libes photo.2image photo.2target photo.camera photo.rectify
photo.target
svn_grass at osgeo.org
svn_grass at osgeo.org
Wed Dec 22 04:44:29 EST 2010
Author: mmetz
Date: 2010-12-22 01:44:29 -0800 (Wed, 22 Dec 2010)
New Revision: 44653
Added:
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_method.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/bilinear.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/bilinear_f.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cubic.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cubic_f.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/nearest.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/readcell.c
Removed:
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_elev.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_files2.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_wind.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/conv.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/matrix.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/perform.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ps_cp.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/rowcol.h
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/write.c
Modified:
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/ls_cameras.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/ls_elev.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/m_mult.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/orthophoto.h
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/orthoref.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2image/group.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/analyze.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/ask.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/equ.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/group.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/zoom.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.camera/main.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_files.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/aver_z.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cp.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/defs.h
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/equ.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/exec.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/get_wind.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/global.h
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/local_proto.h
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/main.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/rectify.c
grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.target/ask_target.c
Log:
merge r44349,44350,44472-44477,44500 from devbr6
Property changes on: grass/branches/releasebranch_6_4/imagery/i.ortho.photo
___________________________________________________________________
Added: svn:mergeinfo
+ /grass/branches/develbranch_6/imagery/i.ortho.photo:44349-44350,44472-44477,44500
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/ls_cameras.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/ls_cameras.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/ls_cameras.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -36,7 +36,7 @@
strcat(buf, ";ls");
if (!full)
strcat(buf, " -C");
- if (ls = popen(buf, "r")) {
+ if ((ls = popen(buf, "r"))) {
while (G_getl(buf, sizeof buf, ls)) {
any = 1;
fprintf(temp, "%s", buf);
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/ls_elev.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/ls_elev.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/ls_elev.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -33,7 +33,7 @@
G__file_name(buf + strlen(buf), element, " ", " ");
strcat(buf, ";ls");
strcat(buf, " -C");
- if (ls = popen(buf, "r")) {
+ if ((ls = popen(buf, "r"))) {
while (G_getl(buf, sizeof buf, ls)) {
any = 1;
fprintf(temp, "%s", buf);
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/m_mult.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/m_mult.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/m_mult.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -11,7 +11,6 @@
{
register int i, j, k, nr, nc, ncols;
char message[256];
- static MATRIX m;
if (a->nrows == 0)
return error("*: arg1 not defined\n");
@@ -31,13 +30,12 @@
nc = b->ncols;
for (i = 0; i < nr; i++)
for (j = 0; j < nc; j++) {
- m.x[i][j] = 0.0;
+ c->x[i][j] = 0.0;
for (k = 0; k < ncols; k++)
- m.x[i][j] += (a->x[i][k] * b->x[k][j]);
+ c->x[i][j] += (a->x[i][k] * b->x[k][j]);
}
- m.nrows = nr;
- m.ncols = nc;
- m_copy(c, &m);
+ c->nrows = nr;
+ c->ncols = nc;
return 1;
}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/orthophoto.h
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/orthophoto.h 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/orthophoto.h 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,5 +1,6 @@
#include <grass/gis.h>
#include <grass/imagery.h>
+#include "mat.h"
/* #define DEBUG 1 */
@@ -57,6 +58,8 @@
int *status;
};
+/* Ortho_Control_Points is identical to Ortho_Photo_Points
+ * Why ? */
struct Ortho_Control_Points
{
int count;
@@ -104,6 +107,7 @@
int con_equation_stat;
double E12[3], N12[3], E21[3], N21[3], Z12[3], Z21[3];
double XC, YC, ZC, omega, phi, kappa;
+ MATRIX M, MI;
};
/* conz_points.c */
@@ -121,13 +125,13 @@
struct Ortho_Camera_File_Ref *,
struct Ortho_Camera_Exp_Init *, double *,
double *, double *, double *, double *,
- double *);
+ double *, MATRIX *, MATRIX *);
int I_ortho_ref(double, double, double, double *, double *, double *,
struct Ortho_Camera_File_Ref *, double, double, double,
- double, double, double);
+ MATRIX);
int I_inverse_ortho_ref(double, double, double, double *, double *, double *,
struct Ortho_Camera_File_Ref *, double, double,
- double, double, double, double);
+ double, MATRIX);
/* ref_points.c */
int I_new_ref_point(struct Ortho_Photo_Points *, double, double, double,
double, int);
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/orthoref.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/orthoref.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/libes/orthoref.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -32,29 +32,26 @@
#include <math.h>
#include "orthophoto.h"
#include <grass/imagery.h>
-#include "mat.h"
#include "matrixdefs.h"
#include "local_proto.h"
-#define MAX_ITERS 10 /* Max iteration is normal equation solution */
-#define CONV_LIMIT 1.0 /* meters */
+#define MAX_ITERS 1000 /* Max iteration is normal equation solution */
+#define CONV_LIMIT 0.1 /* meters */
#ifdef DEBUG
FILE *debug;
char msg[120];
#endif
-static int floating_exception;
-static void catch(int);
-
/* Compute the ortho rectification parameters */
/* XC,YC,ZC, Omega, Phi, Kappa */
int I_compute_ortho_equations(struct Ortho_Control_Points *cpz,
struct Ortho_Camera_File_Ref *cam_info,
struct Ortho_Camera_Exp_Init *init_info,
double *XC, double *YC, double *ZC,
- double *Omega, double *Phi, double *Kappa)
+ double *Omega, double *Phi, double *Kappa,
+ MATRIX *MO, MATRIX *MT)
{
MATRIX delta, epsilon, B, BT, C, E, N, CC, NN, UVW, XYZ, M, WT1;
double meanx, meany;
@@ -74,21 +71,16 @@
Q1 = (double)1.0;
- /*
- * floating_exception = 0;
- * sigfpe = signal (SIGFPE, catch);
- */
-
/* DEBUG */
#ifdef DEBUG
debug = fopen("ortho_compute.rst", "w");
if (debug == NULL) {
- sprintf(msg, "Cant open debug file ortho_analyze.rst\n");
+ sprintf(msg, "Cant open debug file ortho_compute.rst\n");
G_fatal_error(msg);
}
#endif
- /* Need 4 active control points */
+ /* Need at least 4 active control points */
active = 0;
for (i = 0; i < cpz->count; i++) {
if (cpz->status[i] > 0)
@@ -148,7 +140,7 @@
UVW.nrows = 3;
UVW.ncols = 1;
zero(&UVW);
- /* Oreintaiton Matrix M=[3,3] functions of (omega,phi,kappa) */
+ /* Orientation Matrix M=[3,3] functions of (omega,phi,kappa) */
M.nrows = 3;
M.ncols = 3;
zero(&M);
@@ -160,7 +152,7 @@
/******************** Start the solution *****************************/
- /* set Xp, Yp, and CFL form cam_info */
+ /* set Xp, Yp, and CFL from cam_info */
Xp = cam_info->Xp;
Yp = cam_info->Yp;
CFL = cam_info->CFL;
@@ -171,7 +163,7 @@
#endif
/* use initial estimates for XC,YC,ZC,omega,phi,kappa
- * and initial standard deviations if proveded by i.ortho.photo
+ * and initial standard deviations if provided by i.ortho.photo
*
* otherwise set from mean value of all active control points
* - init_info is generated by photo.init (file INIT_EXP)
@@ -194,12 +186,15 @@
WT1.x[4][4] = (Q1 / (init_info->phi_var * init_info->phi_var));
WT1.x[5][5] = (Q1 / (init_info->kappa_var * init_info->kappa_var));
}
- else { /* set from mean values of active control points */
-
-
+ else { /* set from mean values of active control points */
+ double meanX, meanY, meanZ;
+ int nz;
+
/* set intial XC and YC from mean values of control points */
meanx = meany = 0;
- n = 0;
+ meanX = meanY = meanZ = 0;
+ *ZC = 0;
+ n = nz = 0;
first = 1;
for (i = 0; i < cpz->count; i++) {
if (cpz->status[i] <= 0)
@@ -219,20 +214,37 @@
x2 = *(cpz->e1);
Y2 = *(cpz->n2);
y2 = *(cpz->n1);
+ Z1 = *(cpz->z2);
dist_photo =
sqrt(((x1 - x2) * (x1 - x2)) + ((y1 - y2) * (y1 - y2)));
dist_grnd =
sqrt(((X1 - X2) * (X1 - X2)) + ((Y1 - Y2) * (Y1 - Y2)));
+
+ /* set initial ZC from: */
+ /* scale ~= dist_photo/dist_grnd ~= (CFL)/(Z1 - Zc) */
+ /* Zc ~= Z1 + CFL(dist_grnd)/(dist_photo) */
+
+ if (dist_grnd != 0 && dist_photo != 0) {
+ meanZ += Z1 + (CFL * (dist_grnd) / (dist_photo));
+ nz++;
+ }
}
n++;
- meanx += *((cpz->e2)++);
- meany += *((cpz->n2)++);
- ((cpz->e1)++);
- ((cpz->n1)++);
+ meanx += *((cpz->e1)++);
+ meany += *((cpz->n1)++);
+ meanX += *((cpz->e2)++);
+ meanY += *((cpz->n2)++);
+ (cpz->z2)++;
}
- *XC = meanx / n;
- *YC = meany / n;
+ *XC = meanX / n;
+ *YC = meanY / n;
+ *ZC = meanZ / (nz - 1);
+ G_debug(2, "XC: %.2f, YC: %.2f, ZC: %.2f", *XC, *YC, *ZC);
+ meanx = (meanx - x1) / (n - 1);
+ meany = (meany - y1) / (n - 1);
+ meanX = (meanX - X1) / (n - 1);
+ meanY = (meany - Y1) / (n - 1);
/* reset pointers */
for (i = 0; i < cpz->count; i++) {
@@ -242,20 +254,19 @@
(cpz->e2)--;
(cpz->n1)--;
(cpz->n2)--;
+ (cpz->z2)--;
}
- /* set initial ZC from: */
- /* scale ~= dist_photo/dist_grnd ~= (CFL)/(Z1 - Zc) */
- /* Zc ~= Z1 + CFL(dist_grnd)/(dist_photo) */
-
- *ZC = Z1 + (CFL * (dist_grnd) / (dist_photo));
-
/* set initial rotations to zero (radians) */
*Omega = *Phi = 0.0;
/* find an initial kappa */
+ /*
kappa1 = atan2((y2 - y1), (x2 - x1));
kappa2 = atan2((Y2 - Y1), (X2 - X1));
+ */
+ kappa1 = atan2((meany - y1), (meanx - x1));
+ kappa2 = atan2((meanY - Y1), (meanX - X1));
*Kappa = (kappa2 - kappa1);
/* set initial weights */
@@ -362,12 +373,12 @@
fprintf(debug, "\n\t\t\tIn Summation count = %d \n", i);
#endif
- x = *((cpz->e1))++;
- y = *((cpz->n1))++;
- z = *((cpz->z1))++;
- X = *((cpz->e2))++;
- Y = *((cpz->n2))++;
- Z = *((cpz->z2))++;
+ x = *((cpz->e1)++);
+ y = *((cpz->n1)++);
+ z = *((cpz->z1)++);
+ X = *((cpz->e2)++);
+ Y = *((cpz->n2)++);
+ Z = *((cpz->z2)++);
if (cpz->status[i] <= 0)
continue;
@@ -389,7 +400,7 @@
mu = XYZ.x[1][0];
nu = XYZ.x[2][0];
- /* Compute image space coordiantes */
+ /* Compute image space coordinates */
m_mult(&M, &XYZ, &UVW);
/* just an abbreviation */
@@ -397,7 +408,6 @@
V = UVW.x[1][0];
W = UVW.x[2][0];
-
/* Form Partial derivatives of Normal Equations */
xbar = x - Xp;
ybar = y - Yp;
@@ -454,7 +464,7 @@
m_mult(&BT, &E, &C);
/* CC += C */
m_add(&CC, &C, &CC);
- } /* end summation loop over all active control points */
+ } /* end summation loop over all active control points */
/* reset pointers */
for (i = 0; i < cpz->count; i++) {
@@ -479,7 +489,6 @@
CC.x[3][0], CC.x[4][0], CC.x[5][0]);
#endif
-
/* Add weigth matrix of unknowns to NN */
m_add(&NN, &WT1, &NN);
/* Solve for delta */
@@ -503,7 +512,7 @@
delta.x[4][0], delta.x[5][0]);
#endif
- } /* end ITERATION loop */
+ } /* end ITERATION loop */
/* This is the solution */
*XC = epsilon.x[0][0];
@@ -513,6 +522,45 @@
*Phi = epsilon.x[4][0];
*Kappa = epsilon.x[5][0];
+ /* Compute Orientation Matrix M from (Omega, Phi, Kappa); */
+ sw = sin(*Omega);
+ cw = cos(*Omega);
+ sp = sin(*Phi);
+ cp = cos(*Phi);
+ sk = sin(*Kappa);
+ ck = cos(*Kappa);
+
+ MO->nrows = 3;
+ MO->ncols = 3;
+ zero(MO);
+
+ MO->x[0][0] = cp * ck;
+ MO->x[0][1] = cw * sk + (sw * sp * ck);
+ MO->x[0][2] = sw * sk - (cw * sp * ck);
+ MO->x[1][0] = -(cp * sk);
+ MO->x[1][1] = cw * ck - (sw * sp * sk);
+ MO->x[1][2] = sw * ck + (cw * sp * sk);
+ MO->x[2][0] = sp;
+ MO->x[2][1] = -(sw * cp);
+ MO->x[2][2] = cw * cp;
+
+ /* Compute Transposed Orientation Matrix (Omega, Phi, Kappa); */
+
+ MT->nrows = 3;
+ MT->ncols = 3;
+ zero(MT);
+
+ /* M Transposed */
+ MT->x[0][0] = cp * ck;
+ MT->x[1][0] = cw * sk + (sw * sp * ck);
+ MT->x[2][0] = sw * sk - (cw * sp * ck);
+ MT->x[0][1] = -(cp * sk);
+ MT->x[1][1] = cw * ck - (sw * sp * sk);
+ MT->x[2][1] = sw * ck + (cw * sp * sk);
+ MT->x[0][2] = sp;
+ MT->x[1][2] = -(sw * cp);
+ MT->x[2][2] = cw * cp;
+
#ifdef DEBUG
fclose(debug);
#endif
@@ -520,24 +568,16 @@
return (1);
}
-static void catch(int n)
-{
- floating_exception = 1;
- signal(n, catch);
-}
-
/* given ground coordinates (e1,n1,z1) and the solution from above */
/* compute the photo coordinate (e2,n2) position */
int I_ortho_ref(double e1, double n1, double z1,
double *e2, double *n2, double *z2,
struct Ortho_Camera_File_Ref *cam_info,
- double XC, double YC, double ZC,
- double Omega, double Phi, double Kappa)
+ double XC, double YC, double ZC, MATRIX M)
{
- MATRIX UVW, XYZ, M;
+ MATRIX UVW, XYZ;
double U, V, W;
double Xp, Yp, CFL;
- double sw, cw, sp, cp, sk, ck;
/* Initialize and zero the matrices */
/* Object Space Coordinates */
@@ -548,35 +588,13 @@
UVW.nrows = 3;
UVW.ncols = 1;
zero(&UVW);
- /* Orientation Matrix */
- M.nrows = 3;
- M.ncols = 3;
- zero(&M);
/************************ Start the work ******************************/
- /* set Xp, Yp, and CFL form cam_info */
+ /* set Xp, Yp, and CFL from cam_info */
Xp = cam_info->Xp;
Yp = cam_info->Yp;
CFL = cam_info->CFL;
- /* Compute Orientation Matrix M from (Omega, Phi, Kappa); */
- sw = sin(Omega);
- cw = cos(Omega);
- sp = sin(Phi);
- cp = cos(Phi);
- sk = sin(Kappa);
- ck = cos(Kappa);
-
- M.x[0][0] = cp * ck;
- M.x[0][1] = cw * sk + (sw * sp * ck);
- M.x[0][2] = sw * sk - (cw * sp * ck);
- M.x[1][0] = -(cp * sk);
- M.x[1][1] = cw * ck - (sw * sp * sk);
- M.x[1][2] = sw * ck + (cw * sp * sk);
- M.x[2][0] = sp;
- M.x[2][1] = -(sw * cp);
- M.x[2][2] = cw * cp;
-
/* ObjSpace (&XYZ, XC,YC,ZC, X,Y,Z); */
XYZ.x[0][0] = e1 - XC;
XYZ.x[1][0] = n1 - YC;
@@ -601,13 +619,11 @@
int I_inverse_ortho_ref(double e1, double n1, double z1,
double *e2, double *n2, double *z2,
struct Ortho_Camera_File_Ref *cam_info,
- double XC, double YC, double ZC,
- double Omega, double Phi, double Kappa)
+ double XC, double YC, double ZC, MATRIX M)
{
- MATRIX UVW, XYZ, M;
+ MATRIX UVW, XYZ;
double lam, mu, nu;
double Xp, Yp, CFL;
- double sw, cw, sp, cp, sk, ck;
/* Initialize and zero matrices */
/* Object Space Coordinates */
@@ -618,36 +634,13 @@
UVW.nrows = 3;
UVW.ncols = 1;
zero(&UVW);
- /* Orientation Matrix */
- M.nrows = 3;
- M.ncols = 3;
- zero(&M);
/********************** Start the work ********************************/
- /* set Xp, Yp, and CFL form cam_info */
+ /* set Xp, Yp, and CFL from cam_info */
Xp = cam_info->Xp;
Yp = cam_info->Yp;
CFL = cam_info->CFL;
- /* Compute Orientation Matrix (Omega, Phi, Kappa); */
- sw = sin(Omega);
- cw = cos(Omega);
- sp = sin(Phi);
- cp = cos(Phi);
- sk = sin(Kappa);
- ck = cos(Kappa);
-
- /* M Transposed */
- M.x[0][0] = cp * ck;
- M.x[1][0] = cw * sk + (sw * sp * ck);
- M.x[2][0] = sw * sk - (cw * sp * ck);
- M.x[0][1] = -(cp * sk);
- M.x[1][1] = cw * ck - (sw * sp * sk);
- M.x[2][1] = sw * ck + (cw * sp * sk);
- M.x[0][2] = sp;
- M.x[1][2] = -(sw * cp);
- M.x[2][2] = cw * cp;
-
/* ImageSpace */
UVW.x[0][0] = e1 - Xp;
UVW.x[1][0] = n1 - Yp;
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2image/group.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2image/group.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2image/group.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -57,9 +57,9 @@
const int *a = aa, *b = bb;
int n;
- if (n =
+ if ((n =
strcmp(group.group_ref.file[*a].mapset,
- group.group_ref.file[*b].mapset))
+ group.group_ref.file[*b].mapset)) != 0)
return n;
return strcmp(group.group_ref.file[*a].name,
group.group_ref.file[*b].name);
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/analyze.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/analyze.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/analyze.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -410,8 +410,7 @@
temp_points.z2[n],
&e1, &n1, &z1,
&group.camera_ref,
- group.XC, group.YC, group.ZC,
- group.omega, group.phi, group.kappa);
+ group.XC, group.YC, group.ZC, group.MI);
/*****
I_inverse_ortho_ref (group.control_points.e1[n],
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/ask.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/ask.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/ask.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -187,7 +187,7 @@
if (fgets(buf, sizeof buf, fd) == NULL
|| sscanf(buf, "%s %s", name, mapset) != 2)
break;
- if (new_mapset = (strcmp(cur_mapset, mapset) != 0)) {
+ if ((new_mapset = (strcmp(cur_mapset, mapset)) != 0)) {
if (line)
line++;
if (col)
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/equ.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/equ.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/equ.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -38,7 +38,9 @@
&group.ZC,
&group.omega,
&group.phi,
- &group.kappa);
+ &group.kappa,
+ &group.M,
+ &group.MI);
return 0;
}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/group.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/group.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/group.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -55,8 +55,9 @@
const int *a = aa, *b = bb;
int n;
- if (n = strcmp(group.group_ref.file[*a].mapset, group.group_ref.file[*b].
- mapset))
+ if ((n =
+ strcmp(group.group_ref.file[*a].mapset,
+ group.group_ref.file[*b].mapset)) != 0)
return n;
return strcmp(group.group_ref.file[*a].name,
group.group_ref.file[*b].name);
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/zoom.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/zoom.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.2target/zoom.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -43,8 +43,7 @@
for (i = 0; i < 3; i++) {
I_inverse_ortho_ref(spx, spy, spz, trx, try, &trz,
&group.camera_ref,
- group.XC, group.YC, group.ZC,
- group.omega, group.phi, group.kappa);
+ group.XC, group.YC, group.ZC, group.MI);
G_debug(2, "target raster: %.0f %.0f", *trx, *try);
get_z_from_cell2(*try, *trx, &spz);
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.camera/main.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.camera/main.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.camera/main.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -32,6 +32,7 @@
#include <grass/glocale.h>
#include "orthophoto.h"
#include "globals.h"
+
int main(int argc, char *argv[])
{
struct GModule *module;
@@ -77,7 +78,7 @@
strcpy(group, group_opt->answer);
- if( !camera_opt->answer ) {
+ if (!camera_opt->answer) {
/* select the camera to use */
if (!I_ask_camera_any(
_("Enter a camera reference file to be used with this imagery group"),
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_elev.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_elev.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_elev.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,99 +0,0 @@
-/* #include "gis.h" */
-#include <stdlib.h>
-#include <string.h>
-#include <grass/vask.h>
-#include "global.h"
-
-#define ZERO_DATA 0
-
-int ask_elev_data(void)
-{
- char buf1[60];
- char buf2[60];
- char units[10];
- char *tl; /* Target Location */
- char *elev_data; /* Tempoary elevation data layer */
- short ok;
- int no_data_value = ZERO_DATA;
-
- tl = (char *)G_calloc(40, sizeof(char));
- elev_data = (char *)G_calloc(40, sizeof(char));
- tl = G_location();
-
- sprintf(elev_data, "ELEV_DATA");
- sprintf(units, "METERS");
- *buf1 = '\0';
-
- ok = 0;
-
- /** Ask user if modification of elevation data is needed **/
- ok = G_yes("\nModify the data used for elevation ? ", ok);
- if (!ok) {
- select_target_env();
- if ((mapset_elev = G_find_cell(elev_layer, "")) == NULL)
- exit(0);
- select_current_env();
- return (0);
- }
-
- ok = 0;
- while (!ok) {
- /* List options on the screen for the user to answer */
- ok = 1;
- V_clear();
- V_line(1, "Please check the G_mapcalc convention:");
- V_line(3,
- "ELEV_DATA = CELL FILE [MAPSET in LOCATION] [MATH EXPERSION][UNITS]");
-
- V_line(5, "ELEV_DATA : ");
- V_line(6, "CELL FILE : ");
- V_line(7, "MAPSET : ");
- V_line(8, "LOCATION : ");
- V_line(9, "MATH EXPRESSION : ");
- V_line(10, "UNITS : ");
- V_line(12, "NO DATA VALUES = :");
-
-
- /* V_ques ( variable, type, row, col, length) ; */
- V_ques(elev_data, 's', 5, 20, 40);
- V_const(elev_layer, 's', 6, 20, 40);
- V_const(mapset_elev, 's', 7, 20, 40);
- V_const(tl, 's', 8, 20, 40);
- V_ques(buf1, 's', 9, 20, 40);
- V_const(units, 's', 10, 20, 10);
- V_ques(&no_data_value, 'i', 12, 20, 10);
-
- V_intrpt_ok();
- if (!V_call())
- exit(1);
-
- ok = 1;
- sprintf(buf2, "Gmapcalc %s = 'if(%s, %s %s , %d)'", elev_data,
- elev_layer, elev_layer, buf1, no_data_value);
-
- fprintf(stderr,
- "\n\n The following G_mapcalc syntax will be used \n for the modified elevation data\n\n");
- fprintf(stderr, "%s = 'if(%s, %s %s , %d)'", elev_data, elev_layer,
- elev_layer, buf1, no_data_value);
- ok = G_yes("\nDo you accept this G_mapcalc convention \n", ok);
- if (ok) {
-
- /** Set LOCATION to target location **/
- G_setenv("LOCATION_NAME", tl);
-
- /** system GMAPCALC **/
- G_system(buf2);
- /* elev_data becomes the new elevation layer */
- strcpy(elev_layer, elev_data);
- /* need mapset if changed */
- if ((mapset_elev = G_find_cell(elev_layer, "")) == NULL)
- exit(0);
- select_current_env();
-
- /** reset LOCATION to current location **/
- G_setenv("LOCATION_NAME", G_location());
-
- }
- }
- return (0);
-}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_files.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_files.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_files.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -9,15 +9,10 @@
int ask_files(char *groupname)
{
- char line[NFILES][75];
char result[NFILES][15];
- char *err[NFILES];
int repeat;
- int duplicate;
- int list;
- char *name, *mapset;
- int i, k, f1, f2, ln;
- int any;
+ int do_all;
+ int i, k, f1;
int *r;
char **nm;
@@ -25,134 +20,86 @@
nm = new_name;
repeat = 0;
- for (i = 0; i < NFILES; i++)
- err[i] = "";
- f1 = 0;
- for (f2 = f1; f1 < group.group_ref.nfiles; f1 = f2) {
- any = 0;
- ln = 2;
- V_clear();
- V_line(0,
- _("Please select the file(s) you wish to rectify by naming an output file"));
- if (!repeat)
- for (i = 0; i < NFILES; i++)
- result[i][0] = 0;
+
+ /* name extension for rectified maps */
+ extension = (char *)G_malloc(GNAME_MAX * sizeof(char));
+ sprintf(extension, ".ortho");
+
+ repeat = 1;
+ while (repeat) {
repeat = 0;
- for (i = 0; f2 < group.group_ref.nfiles && i < NFILES; i++, f2++) {
- name = group.group_ref.file[f2].name;
- mapset = group.group_ref.file[f2].mapset;
- if (G_find_cell(name, mapset)) {
- sprintf(line[i], "%s in %s", name, mapset);
- dots(line[i], 36);
- V_line(ln, line[i]);
- V_ques(result[i], 's', ln, 37, 14);
- V_const(err[i], 's', ln, 53, 25);
- any = 1;
- ln++;
- }
- }
- if (!any)
- break;
- V_line(ln + 2,
- "(enter list by any name to get a list of existing raster maps)");
+ V_clear();
+ V_line(1, _("Enter an extension to be appended to rectified maps:"));
+ V_ques(extension, 's', 3, 0, 20);
V_intrpt_ok();
if (!V_call())
exit(0);
- /* check the files for illegal names and duplicate names */
-
- list = 0;
- duplicate = 0;
- for (i = 0; i < NFILES; i++) {
- err[i] = "";
- G_strip(result[i]);
- if (result[i][0]) {
- if (strcmp(result[i], "list") == 0) {
- list = 1;
- result[i][0] = 0;
- }
- else if (G_legal_filename(result[i]) < 0) {
- err[i] = "** illegal name **";
- repeat = 1;
- }
- else {
- for (k = 0; k < i; k++) {
- if (!result[k][0])
- continue;
- if (strcmp(result[k], result[i]) != 0)
- continue;
- err[i] = "** duplicate name **";
- duplicate = 1;
- break;
- }
- if (duplicate)
- continue;
- for (k = 0; k < group.group_ref.nfiles; k++) {
- if (ref_list[k] < 0)
- continue;
- if (strcmp(new_name[k], result[i]) != 0)
- continue;
- err[i] = "** duplicate name **";
- duplicate = 1;
- break;
- }
- }
+ /* test for legal file name */
+ sprintf(result[0], "%s%s", group.group_ref.file[0].name, extension);
+ if (G_legal_filename(result[0]) < 0) {
+ G_clear_screen();
+ fprintf(stderr, _("Extension <%s> is illegal"), extension);
+ repeat = G_yes(_("\nChoose another extension? "), 1);
+ if (!repeat) {
+ fprintf(stderr,_("Orthorectification cancelled."));
+ exit(0);
}
}
- if (duplicate)
- repeat = 1;
+ }
+
+ G_debug(1, "Extension: %s", extension);
- /* list the raster maps in the target location. must switch
- * environments to do this
- */
- if (list) {
- repeat = 1;
- select_target_env();
- G_set_list_hit_return(1);
- G_list_element("cell", "raster", G_mapset(), (int (*)())0);
- select_current_env();
- }
+ /* rectify all files ? */
+ do_all = 1;
+ G_clear_screen();
+ do_all = G_yes(_("\nRectify all files in the group? "), do_all);
- if (repeat) {
- f2 = f1;
- continue;
+ /* create list of files to be rectified */
+ f1 = 0;
+ for (i = 0; i < group.group_ref.nfiles && i < NFILES; i++) {
+ int ok = 1;
+ char buf[100];
+
+ if (!do_all) {
+ sprintf(buf, _("\nRectify image <%s>? "), group.group_ref.file[i].name);
+ ok = G_yes(buf, ok);
}
- /* check for existing raster maps
- * this check must occur in the target location, so we switch
- * environments to be in the target location
- */
- select_target_env();
- for (i = 0; i < NFILES; i++) {
- if (result[i][0] && G_find_cell(result[i], G_mapset())) {
- if (!repeat++) {
- repeat = 1;
- fprintf(stderr, "\n");
- fprintf(stderr,
- "** The following raster maps already exist in\n");
- fprintf(stderr, "** LOCATION %s, MAPSET %s:\n\n",
- G_location(), G_mapset());
- }
- fprintf(stderr, "%-18s%s", result[i],
- repeat % 4 ? " " : "\n");
- err[i] = "** file exists **";
- }
+ if (ok) {
+ sprintf(result[i], "%s%s", group.group_ref.file[i].name, extension);
+ *r++ = f1++;
+ *nm++ = G_store(result[i]);
}
- select_current_env();
- if (repeat) {
- repeat = !G_yes("\n\nOk to overwrite? ", 0);
- }
- if (repeat) {
- f2 = f1;
- continue;
- }
+ }
+ for (i = f1; i < NFILES; i++) {
+ result[i][0] = 0;
+ }
- for (i = 0; i < NFILES; i++) {
- if (result[i][0]) {
- *r++ = f1 + i;
- *nm++ = G_store(result[i]);
+ /* check if raster exists in target location/mapset */
+ select_target_env();
+ repeat = 0;
+ G_clear_screen();
+ for (i = 0; i < NFILES; i++) {
+ if (result[i][0] && G_find_cell(result[i], G_mapset())) {
+ if (!repeat++) {
+ repeat = 1;
+ fprintf(stderr, "\n");
+ fprintf(stderr,
+ "** The following raster maps already exist in\n");
+ fprintf(stderr, "** LOCATION %s, MAPSET %s:\n\n",
+ G_location(), G_mapset());
}
+ fprintf(stderr, "%-18s\n", result[i]);
}
}
+ select_current_env();
+ if (repeat) {
+ if (!G_yes("\n\nOk to overwrite? ", 0)) {
+ fprintf(stderr,_("Orthorectification cancelled."));
+ exit(0);
+ }
+ }
+
for (k = 0; k < group.group_ref.nfiles; k++)
if (ref_list[k] >= 0)
return 1;
@@ -160,22 +107,3 @@
G_sleep(3);
exit(0);
}
-
-int dots(char *buf, int n)
-{
- int k;
-
- for (k = 0; *buf; k++)
- buf++;
- if (k >= n)
- return 1;
- *buf++ = ' ';
- k++;
- while (k < n) {
- *buf++ = k % 2 ? '.' : ' ';
- k++;
- }
- *buf = 0;
-
- return 0;
-}
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_files2.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_files2.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_files2.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,80 +0,0 @@
-#include <stdlib.h>
-#include <string.h>
-#include "global.h"
-#include <grass/vask.h>
-
-#define NFILES 15
-
-int ask_file_from_list(char *name, char *mapset)
-{
- char line[NFILES][75];
- char use[NFILES][3];
- int index[NFILES];
- int i, k, n, f1, f2, ln;
- int count;
- int same;
- struct Cell_head win1, win2;
-
- count = 0;
- same = 1;
- for (f1 = 0; f1 < group.group_ref.nfiles; f1++) {
- if (ref_list[f1] >= 0) {
- if (count++ == 0) {
- f2 = ref_list[f1];
- G_get_cellhd(group.group_ref.file[f2].name,
- group.group_ref.file[f2].mapset, &win1);
- }
- else if (same) {
- k = ref_list[f1];
- G_get_cellhd(group.group_ref.file[k].name,
- group.group_ref.file[k].mapset, &win2);
- if (win1.north != win2.north || win1.south != win2.south ||
- win1.east != win2.east || win1.west != win2.west ||
- win1.ns_res != win2.ns_res || win1.ew_res != win2.ew_res)
- same = 0;
- }
- }
- }
- if (count == 1 || same) {
- strcpy(name, group.group_ref.file[f2].name);
- strcpy(mapset, group.group_ref.file[f2].mapset);
- return 1;
- }
- while (1) {
- f1 = 0;
- for (f2 = f1; f1 < group.group_ref.nfiles; f1 = f2) {
- ln = 3;
- V_clear();
- V_line(0,
- "Please mark one file to use as a reference for the window");
- for (i = 0; i < NFILES; i++)
- use[i][0] = 0;
- for (i = 0; f2 < group.group_ref.nfiles && i < NFILES; i++, f2++) {
- if (ref_list[f2] < 0) {
- i--;
- continue;
- }
- n = index[i] = ref_list[f2];
- sprintf(line[i], " %s in %s", group.group_ref.file[n].name,
- group.group_ref.file[n].mapset);
- V_line(ln, line[i]);
- V_ques(use[i], 's', ln, 1, 1);
- ln++;
- }
- if (i == 0)
- break;
- V_intrpt_ok();
- if (!V_call())
- exit(0);
-
- for (i = 0; i < NFILES; i++) {
- if (use[i][0]) {
- n = index[i];
- strcpy(name, group.group_ref.file[n].name);
- strcpy(mapset, group.group_ref.file[n].mapset);
- return 0;
- }
- }
- }
- }
-}
Copied: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_method.c (from rev 44350, grass/branches/develbranch_6/imagery/i.ortho.photo/photo.rectify/ask_method.c)
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_method.c (rev 0)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_method.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -0,0 +1,109 @@
+#include <stdlib.h>
+#include <string.h>
+#include "global.h"
+
+int ask_method(void)
+{
+ int f1, f2, count, max_rows, max_cols;
+ int nx, ny;
+ struct Cell_head win;
+ double max_mb_img, max_mb_elev, max_mb;
+
+ fprintf(stderr, "\n\n");
+ while (1) {
+ char buf[100];
+
+ G_clear_screen();
+ fprintf(stderr, _("Please select one of the following interpolation methods\n"));
+ fprintf(stderr, _(" 1. nearest neighbor\n"));
+ fprintf(stderr, _(" 2. bilinear\n"));
+ fprintf(stderr, _(" 3. bicubic\n"));
+ fprintf(stderr, _(" 4. bilinear with fallback\n"));
+ fprintf(stderr, _(" 5. bicubic with fallback\n"));
+ fprintf(stderr, "> ");
+ if (!G_gets(buf))
+ continue;
+ G_strip(buf);
+
+ if (strcmp(buf, "1") == 0) {
+ interpolate = p_nearest;
+ method = "nearest";
+ break;
+ }
+ if (strcmp(buf, "2") == 0) {
+ interpolate = p_bilinear;
+ method = "bilinear";
+ break;
+ }
+ if (strcmp(buf, "3") == 0) {
+ interpolate = p_cubic;
+ method = "bicubic";
+ break;
+ }
+ if (strcmp(buf, "4") == 0) {
+ interpolate = p_bilinear_f;
+ method = "bilinear_f";
+ break;
+ }
+ if (strcmp(buf, "5") == 0) {
+ interpolate = p_cubic_f;
+ method = "bicubic_f";
+ break;
+ }
+ }
+
+ count = max_rows = max_cols = 0;
+ for (f1 = 0; f1 < group.group_ref.nfiles; f1++) {
+ if (ref_list[f1] >= 0) {
+ f2 = ref_list[f1];
+ G_get_cellhd(group.group_ref.file[f2].name,
+ group.group_ref.file[f2].mapset, &win);
+ if (max_rows < win.rows)
+ max_rows = win.rows;
+ if (max_cols < win.cols)
+ max_cols = win.cols;
+ }
+ }
+
+ ny = (max_rows + BDIM - 1) / BDIM;
+ nx = (max_cols + BDIM - 1) / BDIM;
+
+ max_mb_img = ((double)nx * ny * sizeof(block)) / (1<<20);
+
+ ny = (target_window.rows + BDIM - 1) / BDIM;
+ nx = (target_window.cols + BDIM - 1) / BDIM;
+
+ max_mb_elev = ((double)nx * ny * sizeof(block)) / (1<<20);
+ max_mb = max_mb_img + max_mb_elev + 0.5; /* + 0.5 for rounding */
+ if (max_mb < 1)
+ max_mb = 1;
+
+ fprintf(stderr, "\n\n");
+ while (1) {
+ char buf[100];
+ int seg_mb;
+
+ fprintf(stderr, _("Amount of memory to use in MB\n"));
+ fprintf(stderr, _("RETURN use %d MB to keep all data in RAM\n"), (int)(max_mb));
+ fprintf(stderr, "> ");
+ if (!G_gets(buf))
+ continue;
+
+ if (*buf == 0) { /* all in memory */
+ seg_mb_elev = max_mb_elev;
+ seg_mb_img = max_mb_img;
+ break;
+ }
+
+ G_strip(buf);
+ if ((seg_mb = atoi(buf)) > 0) {
+ seg_mb_elev = seg_mb * max_mb_elev / (max_mb_img + max_mb_elev);
+ seg_mb_img = seg_mb * max_mb_img / (max_mb_img + max_mb_elev);
+ break;
+ }
+ }
+
+ fprintf(stderr, "\n\n");
+
+ return 0;
+}
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_wind.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_wind.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ask_wind.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,205 +0,0 @@
-#include <stdlib.h>
-#include <string.h>
-#include <grass/gis.h>
-#include <grass/vask.h>
-
-static int round(double *);
-static int visually_equal(double, double);
-
-int ask_window(struct Cell_head *window)
-{
- double south, west;
- char buff[64];
- short ok;
- struct Cell_head minimal;
-
-
-
- round(&window->north);
- round(&window->south);
- round(&window->west);
- round(&window->east);
- round(&window->ew_res);
- round(&window->ns_res);
- G_copy(&minimal, window, sizeof(minimal));
-
- window->rows = 0;
- window->cols = 0;
-
- ok = 0;
- while (!ok) {
- /* List window options on the screen for the user to answer */
-
- ok = 1;
- V_clear();
- V_line(0, "Please set the target window");
- V_line(2,
- " ============================= MINIMAL WINDOW ========");
- V_line(3,
- " | North: |");
- V_line(4,
- " | |");
- V_line(5,
- " | ======= GEOREF WINDOW ======= |");
- V_line(6,
- " | | NORTH EDGE: | |");
- V_line(7,
- " | | | |");
- V_line(8,
- " West |WEST EDGE | |EAST EDGE | East");
- V_line(9,
- " | | | |");
- V_line(10,
- " | | SOUTH EDGE: | |");
- V_line(11,
- " | ============================= |");
- V_line(12,
- " | |");
- V_line(13,
- " | South: |");
- V_line(14,
- " =====================================================");
-
- V_line(16,
- " Minimal GRID RESOLUTION Window ");
- V_line(17,
- " --- East-West --- ");
- V_line(18,
- " -- North-South -- ");
- V_line(20,
- "(Minimal window is just large enough to hold entire image)");
-
- /* V_ques ( variable, type, row, col, length) ; */
- V_ques(&window->north, 'd', 6, 36, 11);
- V_ques(&window->south, 'd', 10, 36, 11);
- V_ques(&window->west, 'd', 9, 12, 11);
- V_ques(&window->east, 'd', 9, 52, 11);
- V_ques(&window->ew_res, 'd', 17, 47, 7);
- V_ques(&window->ns_res, 'd', 18, 47, 7);
-
- V_const(&minimal.north, 'd', 3, 36, 11);
- V_const(&minimal.south, 'd', 13, 36, 11);
- V_const(&minimal.west, 'd', 9, 1, 11);
- V_const(&minimal.east, 'd', 9, 66, 11);
- V_const(&minimal.ew_res, 'd', 17, 19, 7);
- V_const(&minimal.ns_res, 'd', 18, 19, 7);
-
- V_intrpt_ok();
- if (!V_call())
- exit(1);
-
- round(&window->north);
- round(&window->south);
- round(&window->east);
- round(&window->west);
- round(&window->ew_res);
- round(&window->ns_res);
-
- if ((window->ns_res <= 0) || (window->ew_res <= 0)) {
- fprintf(stderr, "Illegal resolution value(s)\n");
- ok = 0;
- }
- if (window->north <= window->south) {
- fprintf(stderr, "North must be larger than south\n");
- ok = 0;
- }
- if (window->east <= window->west) {
- fprintf(stderr, "East must be larger than west\n");
- ok = 0;
- }
- if (!ok) {
- fprintf(stderr, "hit RETURN -->");
- G_gets(buff);
- continue;
- }
-
- /* if the north-south is not multiple of the resolution,
- * round the south downward
- */
- south = window->south;
- window->rows =
- (window->north - window->south +
- window->ns_res / 2) / window->ns_res;
- window->south = window->north - window->rows * window->ns_res;
-
- /* do the same for the west */
- west = window->west;
- window->cols =
- (window->east - window->west +
- window->ew_res / 2) / window->ew_res;
- window->west = window->east - window->cols * window->ew_res;
-
- fprintf(stderr, "\n\n");
-
- fprintf(stderr, " north: %12.2f\n", window->north);
- fprintf(stderr, " south: %12.2f", window->south);
-
- if (!visually_equal(window->south, south))
- fprintf(stderr, " (Changed to match resolution)");
- fprintf(stderr, "\n");
-
- fprintf(stderr, " east: %12.2f\n", window->east);
- fprintf(stderr, " west: %12.2f", window->west);
-
- if (!visually_equal(window->west, west))
- fprintf(stderr, " (Changed to match resolution)");
- fprintf(stderr, "\n");
-
- fprintf(stderr, "\n");
- fprintf(stderr, " e-w res: %12.2f\n", window->ew_res);
- fprintf(stderr, " n-s res: %12.2f\n", window->ns_res);
- fprintf(stderr, " total rows: %12d\n", window->rows);
- fprintf(stderr, " total cols: %12d\n", window->cols);
- fprintf(stderr, " total cells: %12d\n", window->rows * window->cols);
- fprintf(stderr, "\n");
-
- ok = 1;
- if (window->north > minimal.north) {
- fprintf(stderr,
- "warning - north falls outside the minimal window\n");
- ok = 0;
- }
- if (window->south < minimal.south) {
- fprintf(stderr,
- "warning - south falls outside the minimal window\n");
- ok = 0;
- }
- if (window->east > minimal.east) {
- fprintf(stderr,
- "warning - east falls outside the minimal window\n");
- ok = 0;
- }
- if (window->west < minimal.west) {
- fprintf(stderr,
- "warning - west falls outside the minimal window\n");
- ok = 0;
- }
-
- ok = G_yes("\nDo you accept this window? ", ok);
- }
-
- return 0;
-}
-
-static int visually_equal(double x, double y)
-{
- char xs[40], ys[40];
-
- if (x == y)
- return 1;
-
- sprintf(xs, "%.2f", x);
- sprintf(ys, "%.2f", y);
-
- return strcmp(xs, ys) == 0;
-}
-
-static int round(double *x)
-{
- char xs[40];
-
- sprintf(xs, "%.2f", *x);
- sscanf(xs, "%lf", x);
-
- return 0;
-}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/aver_z.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/aver_z.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/aver_z.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -3,6 +3,7 @@
int get_aver_elev(struct Ortho_Control_Points *cpz, double *aver_z)
{
double meanz;
+ double *cp = cpz->z2;
int i, n;
/* Need 1 control points */
@@ -18,7 +19,7 @@
continue;
n++;
- meanz += *((cpz->z2)++);
+ meanz += *(cp++);
G_debug(3, "In ortho meanz = %f", meanz);
}
@@ -26,12 +27,5 @@
G_debug(1, "In ortho aver_z = %f", *aver_z);
- /* reset pointers */
- for (i = 0; i < cpz->count; i++) {
- if (cpz->status[i] <= 0)
- continue;
- (cpz->z2)--;
- }
-
return 0;
}
Copied: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/bilinear.c (from rev 44349, grass/branches/develbranch_6/imagery/i.ortho.photo/photo.rectify/bilinear.c)
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/bilinear.c (rev 0)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/bilinear.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -0,0 +1,59 @@
+/*
+ * Name
+ * bilinear.c -- use bilinear interpolation for given row, col
+ *
+ * Description
+ * bilinear interpolation for the given row, column indices.
+ * If the given row or column is outside the bounds of the input map,
+ * the point in the output map is set to NULL.
+ * If any of the 4 surrounding points to be used in the interpolation
+ * is NULL it is filled with is neighbor value
+ */
+
+#include <math.h>
+#include <grass/gis.h>
+#include "global.h"
+
+void p_bilinear(struct cache *ibuffer, /* input buffer */
+ void *obufptr, /* ptr in output buffer */
+ int cell_type, /* raster map type of obufptr */
+ double *row_idx, /* row index */
+ double *col_idx, /* column index */
+ struct Cell_head *cellhd /* information of output map */
+ )
+{
+ int row; /* row indices for interp */
+ int col; /* column indices for interp */
+ int i, j;
+ DCELL t, u; /* intermediate slope */
+ DCELL result; /* result of interpolation */
+ DCELL c[2][2];
+
+ /* cut indices to integer */
+ row = (int)floor(*row_idx - 0.5);
+ col = (int)floor(*col_idx - 0.5);
+
+ /* check for out of bounds - if out of bounds set NULL value and return */
+ if (row < 0 || row + 1 >= cellhd->rows || col < 0 || col + 1 >= cellhd->cols) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+
+ for (i = 0; i < 2; i++)
+ for (j = 0; j < 2; j++) {
+ const DCELL *cellp = CPTR(ibuffer, row + i, col + j);
+ if (G_is_d_null_value(cellp)) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+ c[i][j] = *cellp;
+ }
+
+ /* do the interpolation */
+ t = *col_idx - 0.5 - col;
+ u = *row_idx - 0.5 - row;
+
+ result = G_interp_bilinear(t, u, c[0][0], c[0][1], c[1][0], c[1][1]);
+
+ G_set_raster_value_d(obufptr, result, cell_type);
+}
Copied: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/bilinear_f.c (from rev 44349, grass/branches/develbranch_6/imagery/i.ortho.photo/photo.rectify/bilinear_f.c)
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/bilinear_f.c (rev 0)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/bilinear_f.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -0,0 +1,49 @@
+/*
+ * Name
+ * bilinear_f.c -- use bilinear interpolation with fallback for given row, col
+ *
+ * Description
+ * bilinear interpolation for the given row, column indices.
+ * If the interpolated value (but not the nearest) is null,
+ * it falls back to nearest neighbor.
+ */
+
+#include <grass/gis.h>
+#include <math.h>
+#include "global.h"
+
+void p_bilinear_f(struct cache *ibuffer, /* input buffer */
+ void *obufptr, /* ptr in output buffer */
+ int cell_type, /* raster map type of obufptr */
+ double *row_idx, /* row index */
+ double *col_idx, /* column index */
+ struct Cell_head *cellhd /* cell header of input layer */
+ )
+{
+ /* start nearest neighbor to do some basic tests */
+ int row, col; /* row/col of nearest neighbor */
+ DCELL *cellp, cell;
+
+ /* cut indices to integer */
+ row = (int)floor(*row_idx);
+ col = (int)floor(*col_idx);
+
+ /* check for out of bounds - if out of bounds set NULL value */
+ if (row < 0 || row >= cellhd->rows || col < 0 || col >= cellhd->cols) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+
+ cellp = CPTR(ibuffer, row, col);
+ /* if nearest is null, all the other interps will be null */
+ if (G_is_d_null_value(cellp)) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+ cell = *cellp;
+
+ p_bilinear(ibuffer, obufptr, cell_type, row_idx, col_idx, cellhd);
+ /* fallback to nearest if bilinear is null */
+ if (G_is_d_null_value(obufptr))
+ G_set_raster_value_d(obufptr, cell, cell_type);
+}
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/conv.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/conv.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/conv.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,49 +0,0 @@
-#include "global.h"
-/* conversion routines to convert from view x,y to cell col,row
- * as well as cell col,row to cell east,north
- */
-int view_to_col(View * view, int x)
-{
- return x - view->cell.left;
-}
-
-int view_to_row(View * view, int y)
-{
- return y - view->cell.top;
-}
-
-int col_to_view(View * view, int col)
-{
- return view->cell.left + col;
-}
-
-int row_to_view(View * view, int row)
-{
- return view->cell.top + row;
-}
-
-/* in these next 2 routines, location determines if we are
- * converting from center of the cell (location == .5)
- * top or left edge (location == 0.0)
- * bottom or right edge (location == 1.0)
- */
-
-double row_to_northing(struct Cell_head *cellhd, int row, double location)
-{
- return cellhd->north - (row + location) * cellhd->ns_res;
-}
-
-double col_to_easting(struct Cell_head *cellhd, int col, double location)
-{
- return cellhd->west + (col + location) * cellhd->ew_res;
-}
-
-double northing_to_row(struct Cell_head *cellhd, double north)
-{
- return (cellhd->north - north) / cellhd->ns_res;
-}
-
-double easting_to_col(struct Cell_head *cellhd, double east)
-{
- return (east - cellhd->west) / cellhd->ew_res;
-}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cp.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cp.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cp.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -12,24 +12,23 @@
/* compute photo coordinates of image control points */
/* I_convert_con_points (group.name, &group.control_points, &group.photo_points, group.E12, group.N12); */
-
- sprintf(msg, "Control Z Point file for group [%s] in [%s] \n \n",
+ sprintf(msg, _("Control Z Point file for group [%s] in [%s] \n \n"),
group.name, G_mapset());
- fprintf(stderr, "Computing equations...\n\n");
+ G_verbose_message(_("Computing equations..."));
Compute_ortho_equation();
switch (group.con_equation_stat) {
case -1:
- strcat(msg, "Poorly placed Control Points!\n");
- strcat(msg, "Can not generate the transformation equation.\n");
- strcat(msg, "Run OPTION 7 again!\n");
+ strcat(msg, _("Poorly placed Control Points!\n"));
+ strcat(msg, _("Can not generate the transformation equation.\n"));
+ strcat(msg, _("Run OPTION 7 again!\n"));
break;
case 0:
- strcat(msg, "No active Control Points!\n");
- strcat(msg, "Can not generate the transformation equation.\n");
- strcat(msg, "Run OPTION 7 !\n");
+ strcat(msg, _("No active Control Points!\n"));
+ strcat(msg, _("Can not generate the transformation equation.\n"));
+ strcat(msg, _("Run OPTION 7 !\n"));
break;
default:
return 1;
@@ -46,21 +45,21 @@
if (!I_get_ref_points(group.name, &group.photo_points))
exit(0);
- sprintf(msg, "Reference Point file for group [%s] in [%s] \n \n",
+ sprintf(msg, _("Reference Point file for group [%s] in [%s] \n \n"),
group.name, G_mapset());
Compute_ref_equation();
switch (group.ref_equation_stat) {
case -1:
- strcat(msg, "Poorly placed Reference Points!\n");
- strcat(msg, "Can not generate the transformation equation.\n");
- strcat(msg, "Run OPTION 5 again!\n");
+ strcat(msg, _("Poorly placed Reference Points!\n"));
+ strcat(msg, _("Can not generate the transformation equation.\n"));
+ strcat(msg, _("Run OPTION 5 again!\n"));
break;
case 0:
- strcat(msg, "No active Reference Points!\n");
- strcat(msg, "Can not generate the transformation equation.\n");
- strcat(msg, "Run OPTION 5!\n");
+ strcat(msg, _("No active Reference Points!\n"));
+ strcat(msg, _("Can not generate the transformation equation.\n"));
+ strcat(msg, _("Run OPTION 5!\n"));
break;
default:
E12a = E12[0];
Copied: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cubic.c (from rev 44349, grass/branches/develbranch_6/imagery/i.ortho.photo/photo.rectify/cubic.c)
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cubic.c (rev 0)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cubic.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -0,0 +1,67 @@
+/*
+ * Name
+ * cubic.c -- use cubic convolution interpolation for given row, col
+ *
+ * Description
+ * cubic returns the value in the buffer that is the result of cubic
+ * convolution interpolation for the given row, column indices.
+ * If the given row or column is outside the bounds of the input map,
+ * the corresponding point in the output map is set to NULL.
+ *
+ * If single surrounding points in the interpolation matrix are
+ * NULL they where filled with their neighbor
+ */
+
+#include <grass/gis.h>
+#include <math.h>
+#include "global.h"
+
+void p_cubic(struct cache *ibuffer, /* input buffer */
+ void *obufptr, /* ptr in output buffer */
+ int cell_type, /* raster map type of obufptr */
+ double *row_idx, /* row index (decimal) */
+ double *col_idx, /* column index (decimal) */
+ struct Cell_head *cellhd /* information of output map */
+ )
+{
+ int row; /* row indices for interp */
+ int col; /* column indices for interp */
+ int i, j;
+ DCELL t, u; /* intermediate slope */
+ DCELL result; /* result of interpolation */
+ DCELL val[4]; /* buffer for temporary values */
+ DCELL cell[4][4];
+
+ /* cut indices to integer */
+ row = (int)floor(*row_idx - 0.5);
+ col = (int)floor(*col_idx - 0.5);
+
+ /* check for out of bounds of map - if out of bounds set NULL value */
+ if (row - 1 < 0 || row + 2 >= cellhd->rows ||
+ col - 1 < 0 || col + 2 >= cellhd->cols) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+
+ for (i = 0; i < 4; i++)
+ for (j = 0; j < 4; j++) {
+ const DCELL *cellp = CPTR(ibuffer, row - 1 + i, col - 1 + j);
+ if (G_is_d_null_value(cellp)) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+ cell[i][j] = *cellp;
+ }
+
+ /* do the interpolation */
+ t = *col_idx - 0.5 - col;
+ u = *row_idx - 0.5 - row;
+
+ for (i = 0; i < 4; i++) {
+ val[i] = G_interp_cubic(t, cell[i][0], cell[i][1], cell[i][2], cell[i][3]);
+ }
+
+ result = G_interp_cubic(u, val[0], val[1], val[2], val[3]);
+
+ G_set_raster_value_d(obufptr, result, cell_type);
+}
Copied: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cubic_f.c (from rev 44349, grass/branches/develbranch_6/imagery/i.ortho.photo/photo.rectify/cubic_f.c)
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cubic_f.c (rev 0)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/cubic_f.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -0,0 +1,54 @@
+/*
+ * Name
+ * cubic_f.c -- use cubic interpolation with fallback for given row, col
+ *
+ * Description
+ * cubic interpolation for the given row, column indices.
+ * If the interpolated value (but not the nearest) is null,
+ * it falls back to bilinear. If that interp value is null,
+ * it falls back to nearest neighbor.
+ */
+
+#include <grass/gis.h>
+#include <math.h>
+#include "global.h"
+
+void p_cubic_f(struct cache *ibuffer, /* input buffer */
+ void *obufptr, /* ptr in output buffer */
+ int cell_type, /* raster map type of obufptr */
+ double *row_idx, /* row index */
+ double *col_idx, /* column index */
+ struct Cell_head *cellhd /* cell header of input layer */
+ )
+{
+ /* start nearest neighbor to do some basic tests */
+ int row, col; /* row/col of nearest neighbor */
+ DCELL *cellp, cell;
+
+ /* cut indices to integer */
+ row = (int)floor(*row_idx);
+ col = (int)floor(*col_idx);
+
+ /* check for out of bounds - if out of bounds set NULL value */
+ if (row < 0 || row >= cellhd->rows || col < 0 || col >= cellhd->cols) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+
+ cellp = CPTR(ibuffer, row, col);
+ /* if nearest is null, all the other interps will be null */
+ if (G_is_d_null_value(cellp)) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+ cell = *cellp;
+
+ p_cubic(ibuffer, obufptr, cell_type, row_idx, col_idx, cellhd);
+ /* fallback to bilinear if cubic is null */
+ if (G_is_d_null_value(obufptr)) {
+ p_bilinear(ibuffer, obufptr, cell_type, row_idx, col_idx, cellhd);
+ /* fallback to nearest if bilinear is null */
+ if (G_is_d_null_value(obufptr))
+ G_set_raster_value_d(obufptr, cell, cell_type);
+ }
+}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/defs.h
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/defs.h 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/defs.h 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,65 +1,28 @@
-#include <curses.h>
-#include "orthophoto.h"
-#include <grass/rowio.h>
-#include "rowcol.h"
+/* cache for raster data, code taken from r.proj */
-/* this is a curses structure */
-typedef struct
-{
- int top, left, bottom, right;
-} Window;
+#define L2BDIM 6
+#define BDIM (1<<(L2BDIM))
+#define L2BSIZE (2*(L2BDIM))
+#define BSIZE (1<<(L2BSIZE))
+#define HI(i) ((i)>>(L2BDIM))
+#define LO(i) ((i)&((BDIM)-1))
-/* this is a graphics structure */
-typedef struct
-{
- int top, bottom, left, right;
- int nrows, ncols;
- struct
- {
- int configured;
- struct Cell_head head;
- char name[30];
- char mapset[30];
- int top, bottom, left, right;
- double ew_res, ns_res; /* original map resolution */
- } cell;
-} View;
+typedef DCELL block[BDIM][BDIM]; /* FCELL sufficient ? */
-typedef struct
+struct cache
{
- int type; /* object type */
- int (*handler) (); /* routine to handle the event */
- char *label; /* label to display if MENU or OPTION */
- int binding; /* OPTION bindings */
- int *status; /* MENU,OPTION status */
- int top, bottom, left, right;
-} Objects;
+ int fd;
+ int stride;
+ int nblocks;
+ block **grid;
+ block *blocks;
+ int *refs;
+};
-typedef struct
-{
- double XT, YT, ZT; /* object space */
- int rowT, colT;
- double xt, yt; /* image space */
- int rowt, colt;
-} Tie_Point;
+typedef void (*func) (struct cache *, void *, int, double *, double *, struct Cell_head *);
-typedef struct
-{
- double E12[3], N12[3], E21[3], N21[3];
-} Patch;
+#define BKIDX(c,y,x) ((y) * (c)->stride + (x))
+#define BKPTR(c,y,x) ((c)->grid[BKIDX((c),(y),(x))])
+#define BLOCK(c,y,x) (BKPTR((c),(y),(x)) ? BKPTR((c),(y),(x)) : get_block((c),BKIDX((c),(y),(x))))
+#define CPTR(c,y,x) (&(*BLOCK((c),HI((y)),HI((x))))[LO((y))][LO((x))])
-
-#define MENU_OBJECT 1
-#define OPTION_OBJECT 2
-#define INFO_OBJECT 3
-#define OTHER_OBJECT 4
-
-
-#define MENU(label,handler,status) \
- {MENU_OBJECT,handler,label,0,status,0,0,0,0}
-#define OPTION(label,binding,status) \
- {OPTION_OBJECT,NULL,label,binding,status,0,0,0,0}
-#define INFO(label,status) \
- {INFO_OBJECT,NULL,label,0,status,0,0,0,0}
-#define OTHER(handler,status) \
- {OTHER_OBJECT,handler,NULL,0,status,0,0,0,0}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/equ.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/equ.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/equ.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -39,7 +39,9 @@
&group.ZC,
&group.omega,
&group.phi,
- &group.kappa);
+ &group.kappa,
+ &group.M,
+ &group.MI);
return 0;
}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/exec.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/exec.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/exec.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -23,16 +23,33 @@
struct History hist;
int colr_ok, hist_ok, cats_ok;
long start_time, rectify_time, compress_time;
+ double aver_z;
+ int elevfd;
+ struct cache *ebuffer;
+ G_debug(1, "Open elevation raster: ");
- /* allocate the output cell matrix */
- cell_buf = (CELL **) G_calloc(NROWS, sizeof(void *));
- n = NCOLS * G_raster_size(map_type);
- for (i = 0; i < NROWS; i++) {
- cell_buf[i] = (void *)G_malloc(n);
- G_set_null_value(cell_buf[i], NCOLS, map_type);
+ /* open elevation raster */
+ select_target_env();
+ G_set_window(&target_window);
+ G_debug(1, "target window: rs=%d cs=%d n=%f s=%f w=%f e=%f\n",
+ target_window.rows, target_window.cols, target_window.north,
+ target_window.south, target_window.west, target_window.east);
+
+ elevfd = G_open_cell_old(elev_layer, mapset_elev);
+ if (elevfd < 0) {
+ G_fatal_error(_("Could not open elevation raster"));
+ return 1;
}
+ G_debug(1, "elev layer = %s mapset elev = %s elevfd = %d",
+ elev_layer, mapset_elev, elevfd);
+ ebuffer = readcell(elevfd, seg_mb_elev, 1);
+ G_close_cell(elevfd);
+ /* get an average elevation of the control points */
+ /* this is used only if target cells have no elevation */
+ get_aver_elev(&group.control_points, &aver_z);
+
/* rectify each file */
for (n = 0; n < group.group_ref.nfiles; n++) {
G_debug(2, "I look for files to ortho rectify");
@@ -60,7 +77,7 @@
G_debug(2, "Starting the rectification...");
- if (rectify(name, mapset, result)) {
+ if (rectify(name, mapset, ebuffer, aver_z, result)) {
G_debug(2, "Done. Writing results...");
select_target_env();
if (cats_ok) {
@@ -75,10 +92,7 @@
G_write_history(result, &hist);
select_current_env();
time(&rectify_time);
- if (compress(result))
- time(&compress_time);
- else
- compress_time = rectify_time;
+ compress_time = rectify_time;
report(name, mapset, result, rectify_time - start_time,
compress_time - rectify_time, 1);
}
@@ -87,6 +101,8 @@
report(name, mapset, result, (long)0, (long)0, 0);
}
}
+ close(ebuffer->fd);
+ G_free(ebuffer);
G_done_msg(" ");
return 0;
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/get_wind.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/get_wind.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/get_wind.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,11 +1,62 @@
#include <stdlib.h>
#include <string.h>
+#include <math.h>
#include "global.h"
+int get_ref_window(struct Cell_head *cellhd)
+{
+ int k, f1, f2;
+ int count;
+ struct Cell_head win;
+
+ /* from all the files in the group, get max extends and min resolutions */
+ count = 0;
+ for (f1 = 0; f1 < group.group_ref.nfiles; f1++) {
+ if (ref_list[f1] >= 0) {
+ if (count++ == 0) {
+ f2 = ref_list[f1];
+ G_get_cellhd(group.group_ref.file[f2].name,
+ group.group_ref.file[f2].mapset, cellhd);
+ }
+ else {
+ k = ref_list[f1];
+ G_get_cellhd(group.group_ref.file[k].name,
+ group.group_ref.file[k].mapset, &win);
+ /* extends */
+ if (cellhd->north < win.north)
+ cellhd->north = win.north;
+ if (cellhd->south > win.south)
+ cellhd->south = win.south;
+ if (cellhd->west > win.west)
+ cellhd->west = win.west;
+ if (cellhd->east < win.east)
+ cellhd->east = win.east;
+ /* resolution */
+ if (cellhd->ns_res > win.ns_res)
+ cellhd->ns_res = win.ns_res;
+ if (cellhd->ew_res > win.ew_res)
+ cellhd->ew_res = win.ew_res;
+ }
+ }
+ }
+
+ /* if the north-south is not multiple of the resolution,
+ * round the south downward
+ */
+ cellhd->rows = (cellhd->north - cellhd->south) /cellhd->ns_res + 0.5;
+ cellhd->south = cellhd->north - cellhd->rows * cellhd->ns_res;
+
+ /* do the same for the west */
+ cellhd->cols = (cellhd->east - cellhd->west) / cellhd->ew_res + 0.5;
+ cellhd->west = cellhd->east - cellhd->cols * cellhd->ew_res;
+
+ return 1;
+}
+
int get_target_window(void)
{
- char name[30], mapset[30];
struct Cell_head cellhd;
+ double res;
fprintf(stderr, "\n\n");
while (1) {
@@ -22,33 +73,48 @@
G_strip(buf);
if (strcmp(buf, "1") == 0) {
-
- /**ask_window (&target_window);**/
return 1;
}
if (strcmp(buf, "2") == 0)
break;
}
- ask_file_from_list(name, mapset);
+
+ /* ask for target resolution */
+ while (1) {
+ char buf[100];
- G_debug(1, "ask_file: %s in %s", name, mapset);
+ fprintf(stderr, "Desired target resolution\n");
+ fprintf(stderr,
+ " RETURN determine automatically\n");
+ fprintf(stderr, "> ");
+ if (!G_gets(buf))
+ continue;
- if (G_get_cellhd(name, mapset, &cellhd) < 0)
- exit(EXIT_FAILURE);
+ if (*buf == 0) { /* determine automatically */
+ res = -1;
+ break;
+ }
+ G_strip(buf);
+ if ((res = atof(buf)) <= 0) {
+ fprintf(stderr, "Resolution must be larger than zero!");
+ G_clear_screen();
+ }
+ else
+ break;
+ }
+
+ /* get reference window: max extend, min resolution */
+ get_ref_window(&cellhd);
+
G_debug(1, "current window: n s = %f %f,", cellhd.north,
cellhd.south);
G_debug(1, "current window: w e = %f %f,", cellhd.west,
cellhd.east);
- georef_window(&cellhd, &target_window);
- ask_window(&target_window);
-
-/**
- if(!G_yes("Would you like this window saved as the window in the target location?\n", -1))
- return 0;
-**/
+ georef_window(&cellhd, &target_window, res);
+
select_target_env();
if (G_put_window(&target_window) >= 0)
fprintf(stderr, "Window Saved!\n");
@@ -56,12 +122,14 @@
return 0;
}
-int georef_window(struct Cell_head *w1, struct Cell_head *w2)
+int georef_window(struct Cell_head *w1, struct Cell_head *w2, double res)
{
- double n, e, z1;
+ double n, e, z1, ad;
double n0, e0;
double aver_z;
- double diffew, diffns;
+ struct _corner {
+ double n, e;
+ } nw, ne, se, sw;
/* get an average elevation from the active control points */
get_aver_elev(&group.control_points, &aver_z);
@@ -71,8 +139,7 @@
I_georef(w1->west, w1->north, &e0, &n0, group.E12, group.N12);
I_inverse_ortho_ref(e0, n0, aver_z, &e, &n, &z1, &group.camera_ref,
- group.XC, group.YC, group.ZC, group.omega, group.phi,
- group.kappa);
+ group.XC, group.YC, group.ZC, group.MI);
G_debug(1, "NORTH WEST CORNER");
G_debug(1, "group.E12 = %f %f %f,", group.E12[0], group.E12[1],
@@ -85,11 +152,12 @@
w2->north = w2->south = n;
w2->west = w2->east = e;
+ nw.n = n;
+ nw.e = e;
I_georef(w1->east, w1->north, &e0, &n0, group.E12, group.N12);
I_inverse_ortho_ref(e0, n0, aver_z, &e, &n, &z1, &group.camera_ref,
- group.XC, group.YC, group.ZC, group.omega, group.phi,
- group.kappa);
+ group.XC, group.YC, group.ZC, group.MI);
G_debug(1, "NORTH EAST CORNER");
G_debug(1, "image x = %f y = %f, photo x = %f y = %f", w1->east,
@@ -97,6 +165,8 @@
G_debug(1, "target x = %f y = %f", e, n);
+ ne.n = n;
+ ne.e = e;
if (n > w2->north)
w2->north = n;
if (n < w2->south)
@@ -108,14 +178,15 @@
I_georef(w1->west, w1->south, &e0, &n0, group.E12, group.N12);
I_inverse_ortho_ref(e0, n0, aver_z, &e, &n, &z1, &group.camera_ref,
- group.XC, group.YC, group.ZC, group.omega, group.phi,
- group.kappa);
+ group.XC, group.YC, group.ZC, group.MI);
G_debug(1, "SOUTH WEST CORNER");
G_debug(1, "image x = %f y = %f, photo x = %f y = %f", w1->west,
w1->south, e0, n0);
G_debug(1, "target x = %f y = %f", e, n);
+ sw.n = n;
+ sw.e = e;
if (n > w2->north)
w2->north = n;
if (n < w2->south)
@@ -127,14 +198,15 @@
I_georef(w1->east, w1->south, &e0, &n0, group.E12, N12);
I_inverse_ortho_ref(e0, n0, aver_z, &e, &n, &z1, &group.camera_ref,
- group.XC, group.YC, group.ZC, group.omega, group.phi,
- group.kappa);
+ group.XC, group.YC, group.ZC, group.MI);
G_debug(1, "SOUTH EAST CORNER");
G_debug(1, "image x = %f y = %f, photo x = %f y = %f", w1->east,
w1->south, e0, n0);
G_debug(1, "target x = %f y = %f", e, n);
+ se.n = n;
+ se.e = e;
if (n > w2->north)
w2->north = n;
if (n < w2->south)
@@ -144,20 +216,63 @@
if (e < w2->west)
w2->west = e;
- w2->ns_res = (w2->north - w2->south) / w1->rows;
- w2->ew_res = (w2->east - w2->west) / w1->cols;
+ /* resolution */
+ if (res > 0)
+ w2->ew_res = w2->ns_res = res;
+ else {
+ /* this results in ugly res values, and ns_res != ew_res */
+ /* and is no good for rotation */
+ /*
+ w2->ns_res = (w2->north - w2->south) / w1->rows;
+ w2->ew_res = (w2->east - w2->west) / w1->cols;
+ */
- /* Miori Luca & Mauro Martinelli, ITC-irst 2003: extend region to
- * avoid cut-off of image edges in mountainous terrain:
- * extend target area by (empirically) 15%
- */
- diffew = (w2->east - w2->west);
- diffns = (w2->north - w2->south);
- w2->east = w2->east + 0.15 * diffew;
- w2->west = w2->west - 0.15 * diffew;
- w2->south = w2->south - 0.15 * diffns;
- w2->north = w2->north + 0.15 * diffns;
+ /* alternative: account for rotation and order > 1 */
+ /* N-S extends along western and eastern edge */
+ w2->ns_res = (sqrt((nw.n - sw.n) * (nw.n - sw.n) +
+ (nw.e - sw.e) * (nw.e - sw.e)) +
+ sqrt((ne.n - se.n) * (ne.n - se.n) +
+ (ne.e - se.e) * (ne.e - se.e))) / (2.0 * w1->rows);
+
+ /* E-W extends along northern and southern edge */
+ w2->ew_res = (sqrt((nw.n - ne.n) * (nw.n - ne.n) +
+ (nw.e - ne.e) * (nw.e - ne.e)) +
+ sqrt((sw.n - se.n) * (sw.n - se.n) +
+ (sw.e - se.e) * (sw.e - se.e))) / (2.0 * w1->cols);
+
+ /* make ew_res = ns_res */
+ w2->ns_res = (w2->ns_res + w2->ew_res) / 2.0;
+ w2->ew_res = w2->ns_res;
+
+ /* nice round values */
+ if (w2->ns_res > 1) {
+ if (w2->ns_res < 10) {
+ /* round to first decimal */
+ w2->ns_res = (int)(w2->ns_res * 10 + 0.5) / 10.0;
+ w2->ew_res = w2->ns_res;
+ }
+ else {
+ /* round to integer */
+ w2->ns_res = (int)(w2->ns_res + 0.5);
+ w2->ew_res = w2->ns_res;
+ }
+ }
+ }
+
+ /* adjust extends */
+ ad = w2->north > 0 ? 0.5 : -0.5;
+ w2->north = (int) (ceil(w2->north / w2->ns_res) + ad) * w2->ns_res;
+ ad = w2->south > 0 ? 0.5 : -0.5;
+ w2->south = (int) (floor(w2->south / w2->ns_res) + ad) * w2->ns_res;
+ ad = w2->east > 0 ? 0.5 : -0.5;
+ w2->east = (int) (ceil(w2->east / w2->ew_res) + ad) * w2->ew_res;
+ ad = w2->west > 0 ? 0.5 : -0.5;
+ w2->west = (int) (floor(w2->west / w2->ew_res) + ad) * w2->ew_res;
+
+ w2->rows = (w2->north - w2->south + w2->ns_res / 2.0) / w2->ns_res;
+ w2->cols = (w2->east - w2->west + w2->ew_res / 2.0) / w2->ew_res;
+
G_debug(1, "FINAL");
G_debug(1, "east = %f \n west = %f \n north = %f \n south = %f",
w2->east, w2->west, w2->north, w2->south);
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/global.h
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/global.h 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/global.h 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,52 +1,29 @@
-/* These next defines determine the size of the sub-window that will
- * be held in memory. Larger values will require
- * more memory (but less i/o) If you increase these values, keep in
- * mind that although you think the i/o will decrease, system paging
- * (which goes on behind the scenes) may actual increase the i/o.
- */
-
#include <grass/imagery.h>
#include <grass/ortholib.h>
+#include <grass/glocale.h>
+#include "orthophoto.h"
#include "defs.h"
-/* the large, the worse the results in mountaneous regions!! 128 is MAX!
- * but: the large, the slower - wants a dynamic implementation - TODO
- * possible solution: ratio local elevation range/camera height = 0.003
- */
-
-#define TIE_ROW_DIST 128
-#define TIE_COL_DIST 128
-
-#define NROWS 128
-#define NCOLS 128
-
-/* do not modify past this point */
-
#ifndef GLOBAL
#define GLOBAL extern
#endif
-#define IDX int
-
/* activate debug in Gmakefile */
#ifdef DEBUG3
GLOBAL FILE *Bugsr;
#endif
-GLOBAL ROWCOL row_map[NROWS][NCOLS];
-GLOBAL ROWCOL col_map[NROWS][NCOLS];
-GLOBAL ROWCOL row_min[NROWS];
-GLOBAL ROWCOL row_max[NROWS];
-GLOBAL ROWCOL row_left[NROWS];
-GLOBAL ROWCOL row_right[NROWS];
-GLOBAL IDX row_idx[NROWS];
-GLOBAL int matrix_rows, matrix_cols;
+GLOBAL func interpolate; /* interpolation routine */
+GLOBAL int seg_mb_img, seg_mb_elev;
+GLOBAL char *method;
GLOBAL int temp_fd;
-GLOBAL RASTER_MAP_TYPE map_type;
GLOBAL CELL **cell_buf;
GLOBAL char *temp_name;
+GLOBAL char *extension;
+GLOBAL double target_res;
+
GLOBAL int *ref_list;
GLOBAL char **new_name;
@@ -58,8 +35,6 @@
GLOBAL struct Ortho_Camera_File_Ref cam_info;
GLOBAL struct Cell_head elevhd;
-GLOBAL DCELL *elevbuf;
-GLOBAL int elevfd;
GLOBAL char *elev_layer;
GLOBAL char *mapset_elev;
@@ -72,6 +47,4 @@
GLOBAL struct Cell_head target_window;
-GLOBAL Tie_Point **T_Point;
-
#include "local_proto.h"
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/local_proto.h
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/local_proto.h 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/local_proto.h 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,32 +1,12 @@
-/* ask_elev.c */
-int ask_elev_data(void);
-
/* ask_files.c */
int ask_files(char *);
-int dots(char *, int);
-/* ask_files2.c */
-int ask_file_from_list(char *, char *);
+/* ask_method.c */
+int ask_method(void);
-/* ask_wind.c */
-int ask_window(struct Cell_head *);
-
/* aver_z.c */
int get_aver_elev(struct Ortho_Control_Points *, double *);
-/* compress.c */
-int compress(char *);
-
-/* conv.c */
-int view_to_col(View *, int);
-int view_to_row(View *, int);
-int col_to_view(View *, int);
-int row_to_view(View *, int);
-double row_to_northing(struct Cell_head *, int, double);
-double col_to_easting(struct Cell_head *, int, double);
-double northing_to_row(struct Cell_head *, double);
-double easting_to_col(struct Cell_head *, double);
-
/* cp.c */
int get_conz_points(void);
int get_ref_points(void);
@@ -48,26 +28,34 @@
/* get_wind.c */
int get_target_window(void);
-int georef_window(struct Cell_head *, struct Cell_head *);
+int georef_window(struct Cell_head *, struct Cell_head *, double);
-/* matrix.c */
-int compute_georef_matrix(struct Cell_head *, struct Cell_head *);
-
-/* perform.c */
-int perform_georef(int, void *rast);
-
-/* ps_cp.c */
-int get_psuedo_control_pt(int, int);
-
/* rectify.c */
-int rectify(char *, char *, char *);
+int rectify(char *, char *, struct cache *, double, char *);
+/* readcell.c */
+struct cache *readcell(int, int, int);
+block *get_block(struct cache *, int);
+
/* report.c */
int report(char *, char *, char *, long, long, int);
/* target.c */
int get_target(char *);
-/* write.c */
-int write_map(char *);
-int write_matrix(int, int);
+/* declare resampling methods */
+/* bilinear.c */
+extern void p_bilinear(struct cache *, void *, int, double *, double *,
+ struct Cell_head *);
+/* cubic.c */
+extern void p_cubic(struct cache *, void *, int, double *, double *,
+ struct Cell_head *);
+/* nearest.c */
+extern void p_nearest(struct cache *, void *, int, double *, double *,
+ struct Cell_head *);
+/* bilinear_f.c */
+extern void p_bilinear_f(struct cache *, void *, int, double *, double *,
+ struct Cell_head *);
+/* cubic_f.c */
+extern void p_cubic_f(struct cache *, void *, int, double *, double *,
+ struct Cell_head *);
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/main.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/main.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/main.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -6,11 +6,12 @@
* Markus Neteler <neteler itc.it>,
* Bernhard Reiter <bernhard intevation.de>,
* Glynn Clements <glynn gclements.plus.com>,
- * Hamish Bowman <hamish_b yahoo.com>
+ * Hamish Bowman <hamish_b yahoo.com>,
+ * Markus Metz
*
* PURPOSE: Rectifies an image by using the image to photo coordinate
* transformation matrix
- * COPYRIGHT: (C) 1999-2008 by the GRASS Development Team
+ * COPYRIGHT: (C) 1999-2010 by the GRASS Development Team
*
* This program is free software under the GNU General Public
* License (>=v2). Read the file COPYING that comes with GRASS
@@ -39,7 +40,6 @@
struct GModule *module;
struct Option *group_opt;
-
/* must run in a term window */
G_putenv("GRASS_UI_TERM", "1");
@@ -62,17 +62,18 @@
strcpy(name, group_opt->answer);
- camera = (char *)G_malloc(40 * sizeof(char));
+ camera = (char *)G_malloc(GNAME_MAX * sizeof(char));
elev_layer = (char *)G_malloc(GNAME_MAX * sizeof(char));
mapset_elev = (char *)G_malloc(GMAPSET_MAX * sizeof(char));
- /* get group ref */
+ /* find group */
strcpy(group.name, name);
if (!I_find_group(group.name))
G_fatal_error(_("Group [%s] not found"), group.name);
/* get the group ref */
- I_get_group_ref(group.name, (struct Ref *)&group.group_ref);
+ if (!I_get_group_ref(group.name, (struct Ref *)&group.group_ref))
+ G_fatal_error(_("Could not read REF file for group [%s]"), group.name);
nfiles = group.group_ref.nfiles;
if (nfiles <= 0)
G_fatal_error(_("No files in this group!"));
@@ -88,10 +89,9 @@
/* ask for files to be rectified */
ask_files(group.name);
-
G_debug(1, "Looking for elevation file in group: <%s>", group.name);
- /* get the block elevation layer raster map in target location */
+ /* get the block elevation layer raster map in target location */
if (!I_get_group_elev(group.name, elev_layer, mapset_elev, tl,
math_exp, units, nd))
G_fatal_error(_("No target elevation model selected for group <%s>"),
@@ -104,7 +104,7 @@
G_get_cellhd(elev_layer, mapset_elev, &elevhd);
select_current_env();
- /** look for camera info for this block **/
+ /** look for camera info for this block **/
if (!I_get_group_camera(group.name, camera))
G_fatal_error(_("No camera reference file selected for group <%s>"),
group.name);
@@ -116,7 +116,7 @@
/* get initial camera exposure station, if any */
if (I_find_initial(group.name)) {
if (!I_get_init_info(group.name, &group.camera_exp))
- G_warning(_("Bad format in initial exposusre station file for group <%s>"),
+ G_warning(_("Bad format in initial exposure station file for group <%s>"),
group.name);
}
@@ -130,15 +130,10 @@
select_current_env();
get_target_window();
- /* modify elevation values if needed */
+ /* ask for interpolation method and amount of memory to be used */
+ ask_method();
- /***
- select_target_env();
- ask_elev_data();
- select_current_env();
- ***/
-
- /* go do it */
+ /* go do it */
exec_rectify();
exit(EXIT_SUCCESS);
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/matrix.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/matrix.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/matrix.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,97 +0,0 @@
-
-/**************************************************************
- * compute_georef_matrix (win1, win2)
- *
- */
-#include <stdlib.h>
-#include "global.h"
-
-static int cmp(const void *, const void *);
-
-int compute_georef_matrix(struct Cell_head *win1, struct Cell_head *win2)
-{
- ROWCOL *rmap, *cmap, rr, cc;
- int nrow1, nrow2;
- int ncol1, ncol2;
- double n1, w1, ns_res1, ew_res1;
- double n2, e2, ns_res2, ew_res2;
- double nx, ex;
- double NX, EX;
- int row, col;
- int min, max;
-
- ns_res1 = win1->ns_res;
- ew_res1 = win1->ew_res;
- nrow1 = win1->rows;
- ncol1 = win1->cols;
- n1 = win1->north;
- w1 = win1->west;
-
- ns_res2 = win2->ns_res;
- ew_res2 = win2->ew_res;
- nrow2 = win2->rows;
- ncol2 = win2->cols;
- n2 = win2->north;
- e2 = win2->west;
- matrix_rows = nrow2;
- matrix_cols = ncol2;
-
- /* georef equation is
- * ex = E21a + E21b * col + E21c * row
- * nx = N21a + N21b * col + N21c * row
- *
- * compute common code (for northing) outside east loop
- */
-
- for (n2 = win2->north, row = 0; row < nrow2; row++, n2 -= ns_res2) {
- rmap = row_map[row];
- cmap = col_map[row];
- min = max = -1;
-
-/* G_debug(3, "\n\t got row = \t%d", row); */
-
- /* compute common code */
- EX = E21a + E21c * row;
- NX = N21a + N21c * row;
-
- for (e2 = win2->west, col = 0; col < ncol2; col++, e2 += ew_res2) {
- /* georef e2,n2 */
- ex = EX + E21b * col;
- nx = NX + N21b * col;
-
- rr = (n1 - nx) / ns_res1;
- if (rr < 0 || rr >= nrow1)
- rr = -1;
- else if (min < 0)
- min = max = rr;
- else if (rr < min)
- min = rr;
- else if (rr > max)
- max = rr;
- *rmap++ = rr;
-
- cc = (ex - w1) / ew_res1;
- if (cc < 0 || cc >= ncol1)
- cc = -1;
- *cmap++ = cc;
-
-/* G_debug(4, "\n\tnx = \t%f \tex = \t%f \n\trr = \t%d \tcc = \t%d", nx, ex,rr, cc); */
- }
-
- row_min[row] = min;
- row_max[row] = max;
- row_left[row] = 0;
- row_right[row] = matrix_cols - 1;
- row_idx[row] = row;
- }
- qsort(row_idx, nrow2, sizeof(IDX), cmp);
-
- return 0;
-}
-
-static int cmp(const void *aa, const void *bb)
-{
- const IDX *a = aa, *b = bb;
-
- return (int)(row_min[*a] - row_min[*b]);
-}
Copied: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/nearest.c (from rev 44349, grass/branches/develbranch_6/imagery/i.ortho.photo/photo.rectify/nearest.c)
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/nearest.c (rev 0)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/nearest.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -0,0 +1,40 @@
+/*
+ * nearest.c - returns the nearest neighbor to a given
+ * x,y position
+ */
+
+#include <math.h>
+#include <grass/gis.h>
+#include "global.h"
+
+void p_nearest(struct cache *ibuffer, /* input buffer */
+ void *obufptr, /* ptr in output buffer */
+ int cell_type, /* raster map type of obufptr */
+ double *row_idx, /* row index in input matrix */
+ double *col_idx, /* column index in input matrix */
+ struct Cell_head *cellhd /* cell header of input layer */
+ )
+{
+ int row, col; /* row/col of nearest neighbor */
+ DCELL *cellp;
+
+ /* cut indices to integer and get nearest cell */
+ /* the row_idx, col_idx correction for bilinear/bicubic does not apply here */
+ row = (int)floor(*row_idx);
+ col = (int)floor(*col_idx);
+
+ /* check for out of bounds - if out of bounds set NULL value */
+ if (row < 0 || row >= cellhd->rows || col < 0 || col >= cellhd->cols) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+
+ cellp = CPTR(ibuffer, row, col);
+
+ if (G_is_d_null_value(cellp)) {
+ G_set_null_value(obufptr, 1, cell_type);
+ return;
+ }
+
+ G_set_raster_value_d(obufptr, *cellp, cell_type);
+}
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/perform.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/perform.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/perform.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,105 +0,0 @@
-
-/*======================================================================
- perform.c --
-
- The actual reading and writing of the source imagery
- into the target imagery.
-
-======================================================================*/
-
-#include "global.h"
-
-static ROWCOL *rmap, *cmap, left, right;
-static int do_cell(int, void *, void *);
-
-int rast_size;
-
-int perform_georef(int infd, void *rast)
-{
- int row;
- int curidx;
- int idx;
- int i;
-
- rast_size = G_raster_size(map_type);
- for (row = 0; row < matrix_rows; row++)
- G_set_null_value(cell_buf[row], matrix_cols, map_type);
-
- curidx = 0;
- while (1) {
- /* find first row */
- while (curidx < matrix_rows) {
- idx = row_idx[curidx];
- row = row_min[idx];
- if (row >= 0)
- break;
- curidx++;
- }
-
- /* G_debug(3, " curidx %d", curidx); */
-
- if (curidx >= matrix_rows)
- break;
-
- /* G_debug(3, "read row %d", row); */
-
- if (G_get_raster_row_nomask
- (infd, G_incr_void_ptr(rast, rast_size), row, map_type) < 0)
- return 0;
-
- for (i = curidx; i < matrix_rows; i++) {
- idx = row_idx[i];
- if (row != row_min[idx])
- break;
-
- /* G_debug(4, " process matrix row %d", idx); */
-
- rmap = row_map[idx];
- cmap = col_map[idx];
- left = row_left[idx];
- right = row_right[idx];
- do_cell(row, G_incr_void_ptr(rast, rast_size), cell_buf[idx]);
-
- row_min[idx]++;
- if (row_min[idx] > row_max[idx])
- row_min[idx] = -1;
- row_left[idx] = left;
- row_right[idx] = right;
- }
- }
-
- return 0;
-}
-
-static int do_cell(int row, void *in, void *out)
-{
- int col;
- void *inptr, *outptr;
-
- for (; left <= right; left++) {
- inptr = G_incr_void_ptr(in, cmap[left] * rast_size);
- outptr = G_incr_void_ptr(out, left * rast_size);
- if (rmap[left] < 0)
- continue;
- if (rmap[left] != row)
- break;
- G_raster_cpy(outptr, inptr, 1, map_type);
- }
- for (; left <= right; right--) {
- inptr = G_incr_void_ptr(in, cmap[right] * rast_size);
- outptr = G_incr_void_ptr(out, right * rast_size);
- if (rmap[right] < 0)
- continue;
- if (rmap[right] != row)
- break;
- G_raster_cpy(outptr, inptr, 1, map_type);
- }
- for (col = left; col <= right; col++) {
- inptr = G_incr_void_ptr(in, cmap[col] * rast_size);
- outptr = G_incr_void_ptr(out, col * rast_size);
- if (rmap[col] == row)
- G_raster_cpy(outptr, inptr, 1, map_type);
- }
-
- return 0;
-}
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ps_cp.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ps_cp.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/ps_cp.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,82 +0,0 @@
-#include <stdlib.h>
-#include <string.h>
-#include <grass/glocale.h>
-#include "global.h"
-
-int get_psuedo_control_pt(int tie_row, int tie_col)
-{
- char msg[200];
- struct Ortho_Photo_Points ps_cp;
- int i, j, k;
-
- G_debug(1, "In ps_cp");
-
- /* allocate psuedo struct, max points are four */
- ps_cp.count = 4;
- ps_cp.e1 = (double *)G_malloc(4 * sizeof(double));
- ps_cp.n1 = (double *)G_malloc(4 * sizeof(double));
- ps_cp.e2 = (double *)G_malloc(4 * sizeof(double));
- ps_cp.n2 = (double *)G_malloc(4 * sizeof(double));
- ps_cp.status = (int *)G_malloc(4 * sizeof(int));
- G_debug(1, "ps_cp allocated");
-
- /* pseudo points are four corners taken from T_Points */
- k = 0;
- for (i = 0; i < 2; i++) {
- for (j = 0; j < 2; j++) {
- ps_cp.e1[k] = T_Point[tie_row + i][tie_col + j].xt;
- ps_cp.n1[k] = T_Point[tie_row + i][tie_col + j].yt;
- ps_cp.e2[k] = (j * ((T_Point[tie_row][tie_col + 1].XT -
- T_Point[tie_row][tie_col].XT) /
- target_window.ew_res));
- ps_cp.n2[k] =
- (i *
- ((T_Point[tie_row][tie_col].YT -
- T_Point[tie_row + 1][tie_col].YT) / target_window.ns_res));
- ps_cp.status[k] = 1;
-
- G_debug(2, "\t k = %d\t i = %d\t j = %d", k, i, j);
- G_debug(2, "\t\t e1[k] = %f", ps_cp.e1[k]);
- G_debug(2, "\t\t n1[k] = %f", ps_cp.n1[k]);
- G_debug(2, "\t\t e2[k] = %f", ps_cp.e2[k]);
- G_debug(2, "\t\t n2[k] = %f", ps_cp.n2[k]);
-
- k++;
- }
- }
-
- G_debug(1, "ps_cp initialized");
-
- switch (I_compute_ref_equations(&ps_cp, E12, N12, E21, N21)) {
- case -1:
- G_debug(1, "\tref_equ: case -1");
- strcat(msg, _("Poorly placed psuedo control points. "));
- strcat(msg, _("Cannot generate the transformation equation."));
- break;
- case 0:
- G_debug(1, "\tref_equ: case 0");
- strcat(msg, _("No active psuedo control points"));
- break;
- default:
- G_debug(1, "\tref equ: case good");
-
- E12a = E12[0];
- E12b = E12[1];
- E12c = E12[2];
- N12a = N12[0];
- N12b = N12[1];
- N12c = N12[2];
- E21a = E21[0];
- E21b = E21[1];
- E21c = E21[2];
- N21a = N21[0];
- N21b = N21[1];
- N21c = N21[2];
-
- G_debug(1, "\t\tE21 = %f\t %f\t %f", E21a, E21b, E21c);
- G_debug(1, "\t\tN21 = %f\t %f\t %f", N21a, N21b, N21c);
-
- return 1;
- }
- G_fatal_error(msg);
-}
Copied: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/readcell.c (from rev 44349, grass/branches/develbranch_6/imagery/i.ortho.photo/photo.rectify/readcell.c)
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/readcell.c (rev 0)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/readcell.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -0,0 +1,132 @@
+/*
+ * readcell.c - reads an entire cell layer into a buffer
+ *
+ */
+
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
+#include <sys/types.h>
+#include <sys/stat.h>
+#include <fcntl.h>
+#include <unistd.h>
+#include <grass/gis.h>
+#include <grass/glocale.h>
+#include "global.h"
+
+struct cache *readcell(int fdi, int size, int target_env)
+{
+ DCELL *tmpbuf;
+ struct cache *c;
+ int nrows;
+ int ncols;
+ int row;
+ char *filename;
+ int nx, ny;
+ int nblocks;
+ int i;
+
+ nrows = G_window_rows();
+ ncols = G_window_cols();
+
+ ny = (nrows + BDIM - 1) / BDIM;
+ nx = (ncols + BDIM - 1) / BDIM;
+
+ if (size > 0)
+ nblocks = size * ((1 << 20) / sizeof(block));
+ else
+ nblocks = (nx + ny) * 2; /* guess */
+
+ if (nblocks > nx * ny)
+ nblocks = nx * ny;
+
+ c = G_malloc(sizeof(struct cache));
+ c->stride = nx;
+ c->nblocks = nblocks;
+ c->grid = (block **) G_calloc(nx * ny, sizeof(block *));
+ c->blocks = (block *) G_malloc(nblocks * sizeof(block));
+ c->refs = (int *)G_malloc(nblocks * sizeof(int));
+
+ if (nblocks < nx * ny) {
+ /* Temporary file must be created in output location */
+ select_target_env();
+ filename = G_tempfile();
+ if (!target_env)
+ select_current_env();
+ c->fd = open(filename, O_RDWR | O_CREAT | O_EXCL, 0600);
+ if (c->fd < 0)
+ G_fatal_error(_("Unable to open temporary file"));
+ remove(filename);
+ }
+ else
+ c->fd = -1;
+
+ G_important_message(_("Allocating memory and reading input map..."));
+ G_percent(0, nrows, 5);
+
+ for (i = 0; i < c->nblocks; i++)
+ c->refs[i] = -1;
+
+ tmpbuf = (DCELL *) G_malloc(nx * sizeof(block));
+
+ for (row = 0; row < nrows; row += BDIM) {
+ int x, y;
+
+ for (y = 0; y < BDIM; y++) {
+ G_percent(row + y, nrows, 5);
+
+ if (row + y >= nrows)
+ break;
+
+ G_get_d_raster_row(fdi, &tmpbuf[y * nx * BDIM], row + y);
+ }
+
+ for (x = 0; x < nx; x++)
+ for (y = 0; y < BDIM; y++)
+ if (c->fd >= 0) {
+ if (write
+ (c->fd, &tmpbuf[(y * nx + x) * BDIM],
+ BDIM * sizeof(DCELL)) < 0)
+ G_fatal_error(_("Error writing segment file"));
+ }
+ else
+ memcpy(&c->blocks[BKIDX(c, HI(row), x)][LO(y)][0],
+ &tmpbuf[(y * nx + x) * BDIM],
+ BDIM * sizeof(DCELL));
+ }
+
+ G_free(tmpbuf);
+
+ if (c->fd < 0)
+ for (i = 0; i < c->nblocks; i++) {
+ c->grid[i] = &c->blocks[i];
+ c->refs[i] = i;
+ }
+
+ return c;
+}
+
+block *get_block(struct cache * c, int idx)
+{
+ int replace = rand() % c->nblocks;
+ block *p = &c->blocks[replace];
+ int ref = c->refs[replace];
+ off_t offset = (off_t) idx * sizeof(DCELL) << L2BSIZE;
+
+ if (c->fd < 0)
+ G_fatal_error(_("Internal error: cache miss on fully-cached map"));
+
+ if (ref >= 0)
+ c->grid[ref] = NULL;
+
+ c->grid[idx] = p;
+ c->refs[replace] = idx;
+
+ if (lseek(c->fd, offset, SEEK_SET) < 0)
+ G_fatal_error(_("Error seeking on segment file"));
+
+ if (read(c->fd, p, sizeof(block)) < 0)
+ G_fatal_error(_("Error writing segment file"));
+
+ return p;
+}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/rectify.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/rectify.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/rectify.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -11,230 +11,115 @@
#include "global.h"
#include "local_proto.h"
-int rectify(char *name, char *mapset, char *result)
+int rectify(char *name, char *mapset, struct cache *ebuffer,
+ double aver_z, char *result)
{
- struct Cell_head cellhd, win;
+ struct Cell_head cellhd;
int ncols, nrows;
int row, col;
- int infd;
- void *rast;
- int x_ties, y_ties;
- int tie_row, tie_col;
- int i;
- double n2, e2, z2;
- double nx, ex, zx;
- int r2, c2;
- double row2, col2;
- double aver_z;
+ int infd, outfd;
+ RASTER_MAP_TYPE map_type;
+ int cell_size;
+ void *trast, *tptr;
+ double n1, e1, z1;
+ double nx, ex, nx1, ex1, zx1;
+ double row_idx, col_idx;
+ struct cache *ibuffer;
- G_debug(1, "Open temp elevation file: ");
+ select_current_env();
+ if (G_get_cellhd(name, mapset, &cellhd) < 0)
+ return 0;
- /* open temporary elevation cell layer */
- select_target_env();
+ /* open the file to be rectified
+ * set window to cellhd first to be able to read file exactly
+ */
+ G_set_window(&cellhd);
+ G_debug(2, "cellhd: rs=%d cs=%d n=%f s=%f w=%f e=%f\n",
+ cellhd.rows, cellhd.cols, cellhd.north,
+ cellhd.south, cellhd.west, cellhd.east);
- /**G_set_window (&elevhd);**/
- G_debug(1, "target window: rs=%d cs=%d n=%f s=%f w=%f e=%f\n",
- target_window.rows, target_window.cols, target_window.north,
- target_window.south, target_window.west, target_window.east);
-
- G_set_window(&target_window);
- elevfd = G_open_cell_old(elev_layer, mapset_elev);
-
- /**G_get_cellhd (elev_layer, mapset_elev, &elevhd);**/
- elevbuf = G_allocate_d_raster_buf(); /* enforce DCELL */
-
- /* get an average elevation of the control points */
- /* this is used only if TIE points are outside of the elev_layer boundary */
- get_aver_elev(&group.control_points, &aver_z);
-
-
- if (elevfd < 0) {
- G_debug(1, "CANT OPEN ELEV");
- G_debug(1, "elev layer = %s mapset elev = %s elevfd = %d",
- elev_layer, mapset_elev, elevfd);
+ infd = G_open_cell_old(name, mapset);
+ if (infd < 0) {
return 0;
}
+ map_type = G_get_raster_map_type(infd);
+ cell_size = G_raster_size(map_type);
+ if (strcmp(method, "nearest") != 0) {
+ map_type = DCELL_TYPE;
+ cell_size = G_raster_size(map_type);
+ }
- G_debug(1, "elev layer = %s mapset elev = %s elevfd = %d",
- elev_layer, mapset_elev, elevfd);
+ ibuffer = readcell(infd, seg_mb_img, 0);
- /* alloc Tie_Points */
- y_ties = (int)(target_window.rows / TIE_ROW_DIST) + 2;
- x_ties = (int)(target_window.cols / TIE_COL_DIST) + 2;
+ G_close_cell(infd); /* (pmx) 17 april 2000 */
- G_debug(1, "Number Tie_Points: y_ties %d \tx_ties %d", y_ties,
- x_ties);
+ G_message(_("Rectify <%s@%s> (location <%s>)"),
+ name, mapset, G_location());
+ select_target_env();
+ G_message(_("into <%s@%s> (location <%s>) ..."),
+ result, G_mapset(), G_location());
- T_Point = (Tie_Point **) G_malloc(y_ties * sizeof(Tie_Point *));
- for (i = 0; i < y_ties; i++)
- T_Point[i] = (Tie_Point *) G_malloc(x_ties * sizeof(Tie_Point));
+ G_set_window(&target_window);
+ nrows = target_window.rows;
+ ncols = target_window.cols;
- /* build Tie_Points */
- nrows = 0;
- for (tie_row = 0; tie_row < y_ties; tie_row++) {
- n2 = target_window.north -
- (tie_row * TIE_ROW_DIST * target_window.ns_res) - 1;
- if (n2 <= target_window.south)
- n2 = target_window.south + 1;
- row2 = northing_to_row(&target_window, n2);
- r2 = (int)row2;
+ outfd = G_open_raster_new(result, map_type);
+ trast = G_allocate_raster_buf(map_type);
- if ((G_get_d_raster_row(elevfd, elevbuf, r2)) < 0) {
- G_debug(1, "ERROR reading elevation layer %s fd = %d : row %d",
- elev_layer, elevfd, r2);
- exit(0);
- }
- ncols = 0;
- for (tie_col = 0; tie_col < x_ties; tie_col++) {
- e2 = target_window.west +
- (tie_col * TIE_COL_DIST * target_window.ew_res) + 1;
- if (e2 >= target_window.east)
- e2 = target_window.east - 1;
+ for (row = 0; row < nrows; row++) {
+ n1 = target_window.north - (row + 0.5) * target_window.ns_res;
- G_debug(5, "Tie_Point \t row %d \tcol %d", tie_row, tie_col);
- G_debug(5, "\t east %f\t north %f", e2, n2);
+ G_percent(row, nrows, 2);
- col2 = easting_to_col(&target_window, e2);
- c2 = (int)col2;
+ G_set_null_value(trast, ncols, map_type);
+ tptr = trast;
+ for (col = 0; col < ncols; col++) {
+ DCELL *zp = CPTR(ebuffer, row, col);
- G_debug(5, "\t\t row2 = %f \t col2 = %f", row2, col2);
- G_debug(5, "\t\t r2 = %d \t c2 = %d", r2, c2);
- G_debug(5, "\t\t elevbuf[c2] = %f", (DCELL) elevbuf[c2]);
-
- /* if target TIE point has no elevation, set to aver_z */
- if (G_is_d_null_value(&elevbuf[c2]))
- z2 = aver_z;
+ e1 = target_window.west + (col + 0.5) * target_window.ew_res;
+
+ /* if target cell has no elevation, set to aver_z */
+ if (G_is_d_null_value(zp)) {
+ G_warning(_("No elevation available at row = %d, col = %d"), row, col);
+ z1 = aver_z;
+ }
else
- z2 = (double)elevbuf[c2];
+ z1 = *zp;
- G_debug(5, "\t\t e2 = %f \t n2 = %f \t z2 = %f", e2, n2, z2);
- G_debug(5, "\t\t XC = %f \t YC = %f \t ZC = %f", group.XC,
- group.YC, group.ZC);
- G_debug(5, "\t\t omega = %f \t phi = %f \t kappa = %f",
- group.omega, group.phi, group.kappa);
+ /* target coordinates e1, n1 to photo coordinates ex1, nx1 */
+ I_ortho_ref(e1, n1, z1, &ex1, &nx1, &zx1, &group.camera_ref,
+ group.XC, group.YC, group.ZC, group.M);
- /* ex, nx: photo coordinates */
- I_ortho_ref(e2, n2, z2, &ex, &nx, &zx, &group.camera_ref,
- group.XC, group.YC, group.ZC, group.omega, group.phi,
- group.kappa);
-
G_debug(5, "\t\tAfter ortho ref (photo cords): ex = %f \t nx = %f",
- ex, nx);
+ ex1, nx1);
- /* ex, nx: relative to (row,col) = 0 */
- I_georef(ex, nx, &ex, &nx, group.E21, group.N21);
+ /* photo coordinates ex1, nx1 to image coordinates ex, nx */
+ I_georef(ex1, nx1, &ex, &nx, group.E21, group.N21);
G_debug(5, "\t\tAfter geo ref: ex = %f \t nx = %f", ex, nx);
- T_Point[tie_row][tie_col].XT = e2;
- T_Point[tie_row][tie_col].YT = n2;
- T_Point[tie_row][tie_col].ZT = z2;
- T_Point[tie_row][tie_col].xt = ex;
- T_Point[tie_row][tie_col].yt = nx;
- }
+ /* convert to row/column indices of source raster */
+ row_idx = (cellhd.north - nx) / cellhd.ns_res;
+ col_idx = (ex - cellhd.west) / cellhd.ew_res;
- } /* end build */
+ /* resample data point */
+ interpolate(ibuffer, tptr, map_type, &row_idx, &col_idx, &cellhd);
- /* close elev layer so we can open the file to be rectified */
- select_target_env();
- if (!G_close_cell(elevfd)) {
- G_debug(1, "Can't close the elev file %s [%s in %s]",
- elev_layer, mapset_elev, G_location());
+ tptr = G_incr_void_ptr(tptr, cell_size);
+ }
+ G_put_raster_row(outfd, trast, map_type);
}
+ G_percent(1, 1, 1);
- /* open the result file into target window
- * this open must be first since we change the window later
- * raster maps open for writing are not affected by window changes
- * but those open for reading are
- *
- * also tell open that raster map will have the same format
- * (ie number of bytes per cell) as the file being rectified
- */
- select_current_env();
- if (G_get_cellhd(name, mapset, &cellhd) < 0)
- return 0;
+ G_free(trast);
- select_target_env();
- G_set_window(&target_window);
- G_set_cell_format(cellhd.format);
+ G_close_cell(outfd);
+ close(ibuffer->fd);
+ G_free(ibuffer);
+ if (G_get_cellhd(result, G_mapset(), &cellhd) < 0)
+ return 0;
- select_current_env();
-
- G_copy(&win, &target_window, sizeof(win));
-
- win.west += win.ew_res / 2;
- ncols = target_window.cols;
- col = 0;
-
- for (tie_col = 0; tie_col < (x_ties - 1); tie_col++) {
- G_debug(3, "Patching column %d:", ncols);
-
- if ((win.cols = ncols) > TIE_COL_DIST)
- win.cols = TIE_COL_DIST;
- win.north = target_window.north - win.ns_res / 2;
- nrows = target_window.rows;
- row = 0;
-
- for (tie_row = 0; tie_row < (y_ties - 1); tie_row++) {
- G_debug(5, "Patching %d row:", nrows);
-
- if ((win.rows = nrows) > TIE_ROW_DIST)
- win.rows = TIE_ROW_DIST;
-
- get_psuedo_control_pt(tie_row, tie_col);
-
- G_debug(5, "\t got psuedo pts: row %d \t col %d", tie_row, tie_col);
-
- compute_georef_matrix(&cellhd, &win);
-
- G_debug(5, "\t\tcompute geo matrix");
-
- /* open the source imagery file to be rectified */
- /* set window to cellhd first to be able to read file exactly */
- select_current_env();
- G_set_window(&cellhd);
- infd = G_open_cell_old(name, mapset);
- if (infd < 0) {
- close(infd);
- return 0;
- }
- map_type = G_get_raster_map_type(infd);
- rast =
- (void *)G_calloc(G_window_cols() + 1,
- G_raster_size(map_type));
- G_set_null_value(rast, G_window_cols() + 1, map_type);
-
- /* perform the actual data rectification */
- perform_georef(infd, rast);
-
- G_debug(5, "\t\tperform georef");
- G_debug(5, "\t\twrite matrix");
-
- /* close the source imagery file and free the buffer */
- select_current_env();
- G_close_cell(infd);
- G_free(rast);
- /* select_current_env(); */
- select_target_env();
-
-
- /* write of the data rectified into the result file */
- write_matrix(row, col);
-
- nrows -= win.rows;
- row += win.rows;
- win.north -= (win.ns_res * win.rows);
- }
-
- ncols -= win.cols;
- col += win.cols;
- win.west += (win.ew_res * win.cols);
- G_percent(col, col + ncols, 1);
- }
-
- select_target_env();
-
if (cellhd.proj == 0) { /* x,y imagery */
cellhd.proj = target_window.proj;
cellhd.zone = target_window.zone;
@@ -252,9 +137,6 @@
name, mapset);
}
- target_window.compressed = cellhd.compressed;
- G_close_cell(infd);
- write_map(result);
select_current_env();
return 1;
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/rowcol.h
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/rowcol.h 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/rowcol.h 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1 +0,0 @@
-#define ROWCOL short
Deleted: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/write.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/write.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.rectify/write.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -1,76 +0,0 @@
-#include <sys/types.h>
-#include <sys/stat.h>
-#include <fcntl.h>
-#include <unistd.h>
-#include <string.h>
-#include <errno.h>
-#include <grass/glocale.h>
-#include "global.h"
-
-
-int write_matrix(int row, int col)
-{
- int n;
-
- select_target_env();
-
- /* create temp file if it doesn't eexist */
- if (!temp_fd || (fcntl(temp_fd, F_GETFD) == -1)) {
- temp_name = G_tempfile();
- temp_fd = creat(temp_name, 0660);
- }
-
- for (n = 0; n < matrix_rows; n++) {
- off_t offset;
-
- offset =
- ((off_t) row++ * target_window.cols +
- col) * G_raster_size(map_type);
- lseek(temp_fd, offset, SEEK_SET);
-
- if (write(temp_fd, cell_buf[n], G_raster_size(map_type) * matrix_cols)
- != G_raster_size(map_type) * matrix_cols) {
- unlink(temp_name);
- G_fatal_error(_("Unable to write temp file: %s"),
- strerror(errno));
- }
- }
-
- select_current_env();
-
- return 0;
-}
-
-
-int write_map(char *name)
-{
- int fd, row;
- void *rast;
-
- G_set_window(&target_window);
-
- rast = G_allocate_raster_buf(map_type);
- close(temp_fd);
- temp_fd = open(temp_name, F_DUPFD);
- fd = G_open_raster_new(name, map_type);
-
- if (fd <= 0)
- G_fatal_error(_("Unable to open map %s"), name);
-
- for (row = 0; row < target_window.rows; row++) {
- if (read(temp_fd, rast, target_window.cols * G_raster_size(map_type))
- != target_window.cols * G_raster_size(map_type))
- G_fatal_error(_("Unable to write row %d"), row);
-
- if (G_put_raster_row(fd, rast, map_type) < 0) {
- G_fatal_error(_("Unable to write raster map. You might want to check available disk space and write permissions."));
- unlink(temp_name);
- }
- }
-
- close(temp_fd);
- unlink(temp_name);
- G_close_cell(fd);
-
- return 0;
-}
Modified: grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.target/ask_target.c
===================================================================
--- grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.target/ask_target.c 2010-12-22 09:38:28 UTC (rev 44652)
+++ grass/branches/releasebranch_6_4/imagery/i.ortho.photo/photo.target/ask_target.c 2010-12-22 09:44:29 UTC (rev 44653)
@@ -110,7 +110,7 @@
fprintf(stderr, "\n");
tot_len = len;
}
- if (ok = (G__mapset_permissions(buf) == 1))
+ if ((ok = (G__mapset_permissions(buf)) == 1))
any_ok = 1;
fprintf(stderr, "%s%-*s", ok ? "(+)" : " ", len, buf);
}
More information about the grass-commit
mailing list