/* Dump out the ascii format DEM file */
-void rawDumpAsciiDEM( fgRAWDEM *raw ) {
+void rawDumpAsciiDEM( fgRAWDEM *raw, char *path, int ilon, int ilat ) {
char outfile[256];
+ char tmp[256];
+ int lon, lat;
+ char lon_sign, lat_sign;
+ int i, j;
+ FILE *fd;
+
/* Generate output file name */
-
-
- /* Dump the "A" record */
- /* get the name field (144 characters) */
- for ( i = 0; i < 144; i++ ) {
- name[i] = fgetc(fd);
+ if ( ilon >= 0 ) {
+ lon = ilon;
+ lon_sign = 'e';
+ } else {
+ lon = -ilon;
+ lon_sign = 'w';
}
- name[i+1] = '\0';
-
- /* clean off the whitespace at the end */
- for ( i = strlen(name)-2; i > 0; i-- ) {
- if ( !isspace(name[i]) ) {
- i=0;
- } else {
- name[i] = '\0';
- }
+
+ if ( ilat >= 0 ) {
+ lat = ilat;
+ lat_sign = 'n';
+ } else {
+ lat = -ilat;
+ lat_sign = 's';
}
- printf(" Quad name field: %s\n", name);
- /* get quadrangle id (now part of previous section */
- /* next_token(fd, dem_quadrangle); */
- /* printf(" Quadrangle = %s\n", dem_quadrangle); */
+ sprintf(outfile, "%s/%c%03d%c%03d.dem", path, lon_sign, lon, lat_sign, lat);
+
+ printf("outfile = %s\n", outfile);
+
+ if ( (fd = fopen(outfile, "w")) == NULL ) {
+ printf("Error opening output file = %s\n", outfile);
+ exit(-1);
+ }
+
+ /* Dump the "A" record */
+
+ /* print descriptive header (144 characters) */
+ sprintf(tmp, "%s - Generated from a 30 arcsec binary DEM", outfile);
+ fprintf(fd, "%-144s", tmp);
/* DEM level code, 3 reflects processing by DMA */
- inum = next_int(fd);
- printf(" DEM level code = %d\n", inum);
+ fprintf(fd, "%6d", 1);
/* Pattern code, 1 indicates a regular elevation pattern */
- inum = next_int(fd);
- printf(" Pattern code = %d\n", inum);
+ fprintf(fd, "%6d", 1);
/* Planimetric reference system code, 0 indicates geographic
* coordinate system. */
- inum = next_int(fd);
- printf(" Planimetric reference code = %d\n", inum);
+ fprintf(fd, "%6d", 0);
/* Zone code */
- inum = next_int(fd);
- printf(" Zone code = %d\n", inum);
+ fprintf(fd, "%6d", 0);
/* Map projection parameters (ignored) */
for ( i = 0; i < 15; i++ ) {
- dnum = next_double(fd);
- /* printf("%d: %f\n",i,dnum); */
+ fprintf(fd, "%6.1f%18s", 0.0, "");
}
- /* Units code, 3 represents arc-seconds as the unit of measure for
+ /* Units code, 3 represents arc-seconds as the unit of measure for
* ground planimetric coordinates throughout the file. */
- inum = next_int(fd);
- if ( inum != 3 ) {
- printf(" Unknown (X,Y) units code = %d!\n", inum);
- exit(-1);
- }
+ fprintf(fd, "%6d", 3);
/* Units code; 2 represents meters as the unit of measure for
* elevation coordinates throughout the file. */
- inum = next_int(fd);
- if ( inum != 2 ) {
- printf(" Unknown (Z) units code = %d!\n", inum);
- exit(-1);
- }
+ fprintf(fd, "%6d", 2);
/* Number (n) of sides in the polygon which defines the coverage of
* the DEM file (usually equal to 4). */
- inum = next_int(fd);
- if ( inum != 4 ) {
- printf(" Unknown polygon dimension = %d!\n", inum);
- exit(-1);
- }
+ fprintf(fd, "%6d", 4);
/* Ground coordinates of bounding box in arc-seconds */
- dem_x1 = m->originx = next_exp(fd);
- dem_y1 = m->originy = next_exp(fd);
- printf(" Origin = (%.2f,%.2f)\n", m->originx, m->originy);
+ fprintf(fd, "%20.15fD+06", ilon * 3600.0 / 1000000.0);
+ fprintf(fd, "%20.15fD+06", ilat * 3600.0 / 1000000.0);
- dem_x2 = next_exp(fd);
- dem_y2 = next_exp(fd);
+ fprintf(fd, "%20.15fD+06", ilon * 3600.0 / 1000000.0);
+ fprintf(fd, "%20.15fD+06", (ilat+1) * 3600.0 / 1000000.0);
- dem_x3 = next_exp(fd);
- dem_y3 = next_exp(fd);
+ fprintf(fd, "%20.15fD+06", (ilon+1) * 3600.0 / 1000000.0);
+ fprintf(fd, "%20.15fD+06", (ilat+1) * 3600.0 / 1000000.0);
- dem_x4 = next_exp(fd);
- dem_y4 = next_exp(fd);
+ fprintf(fd, "%20.15fD+06", (ilon+1) * 3600.0 / 1000000.0);
+ fprintf(fd, "%20.15fD+06", (ilat) * 3600.0 / 1000000.0);
/* Minimum/maximum elevations in meters */
- dem_z1 = next_exp(fd);
- dem_z2 = next_exp(fd);
- printf(" Elevation range %.4f %.4f\n", dem_z1, dem_z2);
+ fprintf(fd, " %20.15E", (double)raw->tmp_min);
+ fprintf(fd, " %20.15E", (double)raw->tmp_max);
/* Counterclockwise angle from the primary axis of ground
* planimetric referenced to the primary axis of the DEM local
* reference system. */
- next_token(fd, token);
+ fprintf(fd, "%6.1f", 0.0);
/* Accuracy code; 0 indicates that a record of accuracy does not
* exist and that no record type C will follow. */
- /* DEM spacial resolution. Usually (3,3,1) (3,6,1) or (3,9,1)
+ fprintf(fd, "%24d", 0);
+
+ /* DEM spacial resolution. Usually (3,3) (3,6) or (3,9)
* depending on latitude */
- /* I will eventually have to do something with this for data at
- * higher latitudes */
- next_token(fd, token);
- printf(" accuracy & spacial resolution string = %s\n", token);
- i = strlen(token);
- printf(" length = %d\n", i);
-
- ptr = token + i - 12;
- printf(" last field = %s = %.2f\n", ptr, atof(ptr));
- ptr[0] = '\0';
-
- ptr = ptr - 12;
- m->col_step = atof(ptr);
- printf(" last field = %s = %.2f\n", ptr, m->row_step);
- ptr[0] = '\0';
-
- ptr = ptr - 12;
- m->row_step = atof(ptr);
- printf(" last field = %s = %.2f\n", ptr, m->col_step);
- ptr[0] = '\0';
-
- /* accuracy code = atod(token) */
- inum = atoi(token);
- printf(" Accuracy code = %d\n", inum);
-
- printf(" column step = %.2f row step = %.2f\n",
- m->col_step, m->row_step);
+ fprintf(fd, "%12.6E", 30.0);
+ fprintf(fd, "%12.6E", 30.0);
+
+ /* accuracy code */
+ fprintf(fd, "%12.6E", 1.0);
+
/* dimension of arrays to follow (1)*/
- next_token(fd, token);
+ fprintf(fd, "%6d", 1);
/* number of profiles */
- dem_num_profiles = m->cols = m->rows = next_int(fd);
- printf(" Expecting %d profiles\n", dem_num_profiles);
+ fprintf(fd, "%6d", 3600 / raw->ydim + 1);
+ /* pad the end */
+ fprintf(fd, "%160s", "");
+
+
+ /* Dump "B" records */
+
+ for ( j = 0; j <= 120; j++ ) {
+ /* row / column id of this profile */
+ fprintf(fd, "%6d%6d", 1, j + 1);
+
+ /* Number of columns and rows (elevations) in this profile */
+ fprintf(fd, "%6d%6d", 3600 / raw->xdim + 1, 1);
+
+ /* Ground planimetric coordinates (arc-seconds) of the first
+ * elevation in the profile */
+ fprintf(fd, "%20.15fD+06", ilon * 3600.0 / 1000000.0);
+ fprintf(fd, "%20.15fD+06", (ilat * 3600.0 + j * raw->ydim) / 1000000.0);
+
+ /* Elevation of local datum for the profile. Always zero for
+ * 1-degree DEM, the reference is mean sea level. */
+ fprintf(fd, "%6.1f", 0.0);
+ fprintf(fd, "%18s", "");
+
+ /* Minimum and maximum elevations for the profile. */
+ fprintf(fd, " %20.15E", 0.0);
+ fprintf(fd, " %20.15E", 0.0);
+
+ /* One (usually) dimensional array (prof_num_cols,1) of elevations */
+ for ( i = 0; i <= 120; i++ ) {
+ fprintf(fd, "%6.0f", raw->edge[j][i]);
+ }
+ }
+
+ fprintf(fd, "\n");
+
+ fclose(fd);
}
/* Read a horizontal strip of (1 vertical degree) from the raw DEM
* file specified by the upper latitude of the stripe specified in
- * degrees */
-void rawReadStrip( fgRAWDEM *raw, int lat_degrees ) {
+ * degrees. The output the individual ASCII format DEM tiles. */
+void rawProcessStrip( fgRAWDEM *raw, int lat_degrees, char *path ) {
int lat, yrange;
int i, j, index, row, col;
- int min = 0, max = 0;
+ int min, max;
int span, num_degrees, tile_width;
int xstart, xend;
/* map ocean to 0 for now */
raw->strip[index][j] = 0;
}
-
- if ( raw->strip[index][j] < min) {
- min = raw->strip[index][j];
- }
-
- if ( raw->strip[index][j] > max) {
- max = raw->strip[index][j];
- }
}
}
- printf("min = %d max = %d\n", min, max);
/* extract individual tiles from the strip */
span = raw->ncols * raw->xdim;
xstart = i * tile_width;
xend = xstart + 120;
+ min = 10000; max = -10000;
for ( row = 0; row < yrange; row++ ) {
for ( col = xstart; col < xend; col++ ) {
/* Copy from strip to pixel centered tile. Yep,
* row/col are reversed here. raw->strip is backwards
* for convenience. I am converting to [x,y] now. */
raw->center[col-xstart][row] = raw->strip[row][col];
+
+ if ( raw->strip[row][col] < min) {
+ min = raw->strip[row][col];
+ }
+
+ if ( raw->strip[row][col] > max) {
+ max = raw->strip[row][col];
+ }
}
}
+ raw->tmp_min = min;
+ raw->tmp_max = max;
+
/* Convert from pixel centered to pixel edge values */
rawConvertCenter2Edge(raw);
/* Dump out the ascii format DEM file */
- rawDumpAsciiDEM(raw);
+ rawDumpAsciiDEM(raw, path, (raw->rootx / 3600) + i, lat_degrees - 1);
}
}
/* $Log$
-/* Revision 1.2 1998/03/03 02:04:01 curt
-/* Starting DEM Ascii format output routine.
+/* Revision 1.3 1998/03/03 13:10:29 curt
+/* Close to a working version.
/*
+ * Revision 1.2 1998/03/03 02:04:01 curt
+ * Starting DEM Ascii format output routine.
+ *
* Revision 1.1 1998/03/02 23:31:01 curt
* Initial revision.
*
projects / flightgear.git / commitdiff