david pratt (pratt++at++medusa.cs.nps.navy.mil)
Wed, 22 Mar 1995 20:47:08 -0800
Dave
--Dave Pratt pratt++at++cs.nps.navy.mil (408) 656-2865 fax (408) 656-2814 Department of Computer Science, Naval Postgraduate School, Monterey, CA 93943 These are my opinions, talk to the PAO for the Navy's.
/* This is the testbed for texturing the terrain in Load_DMA.c, the backup program for non-textured Load_DMA.c is Load_DMA.bak */
#include <stdio.h> #include <pf.h> #define _DMA_LOADER_ #include "dma_dted.h" #include "math.h" // #include "constants.h" #define NUMCOLORS 128
int max_row,max_col; double delta_lng, delta_lat; int origin[86][2]; static int quadr[4][2]; pfVec3 *dirt_grid,*normal_array; void *arena;
pfVec2 texcoords[]={ {0.0f, 0.0f}, {1.0f, 0.0f}, {0.0f, 1.0f}, {1.0f, 1.0f}}; ushort texlist[] = { 0, 1, 2, 1, 3, 2, 0, 1, 2, 1, 3, 2,
0, 1, 2, 1, 3, 2, 0, 1, 2, 1, 3, 2,
0, 1, 2, 1, 3, 2, 0, 1, 2, 1, 3, 2,
0, 1, 2, 1, 3, 2, 0, 1, 2, 1, 3, 2,
0, 1, 2, 1, 3, 2, 0, 1, 2, 1, 3, 2,
0, 1, 2, 1, 3, 2, 0, 1, 2, 1, 3, 2,
0, 1, 2, 1 }; /* 76 elevation points per t-strip, 5776 per gset */
pfGroup *LoadDTED(char *FileName) { int numvert,ix,num_of_points, num_child, tot_child=85, nlev=1, totlev = 5, nquad, nparent; float minelev,maxelev; pfVec4 *color_array; ushort *index_array,*color_index_array; DTED *dted; pfGeoSet *map; int row_cnt,coll_cnt; long *prim_cnt_arr; pfGeode *geode; pfGeoState *gstate; pfGroup *group[86]; pfTexture *tex; pfTexEnv *tev; /* if (argc < 2){ printf("Usage dma_loader <filename>\n"); exit(0); } */
/* get IRIS Performer shared memory arena */ arena = pfGetSharedArena();
gstate = pfNewGState(arena);
/* Set up textures & gstates structures */ tev = pfNewTEnv (arena); tex = pfNewTex (arena); /* define geostate for terrain */ /* if ( pfLoadTexFile (tex, "lig_dirt.rgb")) { pfGStateMode(gstate, PFSTATE_ENLIGHTING,1); printf("Have loaded the texture file\n"); pfGStateMode(gstate, PFSTATE_ENTEXTURE, 1); pfGStateAttr (gstate, PFSTATE_TEXENV, tev); pfGStateAttr (gstate, PFSTATE_TEXTURE, tex); } else { printf("Did NOT Load a texture file\n"); pfGStateMode(gstate, PFSTATE_ENLIGHTING, 0); pfGStateMode(gstate, PFSTATE_ENTEXTURE, 0); pfGStateMode(gstate, PFSTATE_ENWIREFRAME, 0); } */
pfGStateMode(gstate, PFSTATE_ENLIGHTING, 0); pfGStateMode(gstate, PFSTATE_ENTEXTURE, 0); pfGStateMode(gstate, PFSTATE_ENWIREFRAME, 0);
quadr[0][0] = quadr[0][1] = quadr[1][1] = quadr[3][0] = 0; quadr[1][0] = quadr[3][1] = quadr[2][0] = quadr[2][1] = 1;
dted = DMA_read_DTED_terrain(FileName,NULL); min_max_elev(dted, &minelev, &maxelev); /* printf("\n minelev = %f, maxelev = %f\n", minelev, maxelev); */
/* Set smallest quadrangle size to be 76 x 76 data points */ (dted->Lat_Count)[0] = (dted->Lng_Count)[0] = '0'; (dted->Lat_Count)[1] = (dted->Lng_Count)[1] = '0'; (dted->Lat_Count)[2] = (dted->Lng_Count)[2] = '7'; (dted->Lat_Count)[3] = (dted->Lng_Count)[3] = '6';
sscanf(dted->Lat_Count,"%4d",&row_cnt); sscanf(dted->Lng_Count,"%4d",&coll_cnt);
/*build the color ramp*/ color_array = DMA_build_color_ramp(NUMCOLORS);
/*The Normals are next*/ normal_array = DMA_build_normal_array(dted);
for(num_child = 1; num_child <= tot_child; num_child++) { group[num_child] = pfNewGroup();
/* printf("New group %d created\n",num_child); */ /* Compute indices of quadrangle origin */ /* if(NULL == (origin = (int *)pfMalloc(sizeof(int)*2*tot_child,arena))){ printf("Unable to pfMalloc the origin array\n"); exit(0); } */
if(num_child>1) { nparent = (num_child+2)/4; if(((nparent == 2) || (nparent == 6) || (nparent == 22) ) && (((num_child+2)%4) == 0) ) nlev ++; nquad = ((num_child+2)%4); origin[num_child][0] = origin[nparent][0] + quadr[nquad][0]*(int)pow((double)2.0,(double)(totlev-nlev))/2*(coll_cnt-1); origin[num_child][1] = origin[nparent][1] + quadr[nquad][1]*(int)pow((double)2.0,(double)(totlev-nlev))/2*(row_cnt-1); /* printf("nlev = %d, nquad = %d: num_child= %d, xorigin = %d, yorigin = %d\n", nlev, nquad,num_child, origin[num_child][0],origin[num_child][1]); */
} else { origin[num_child][0] = 0 ; origin[num_child][1] = 0 ; }
if(num_child > 21) { geode = pfNewGeode(); for(nquad = 0; nquad < 4; nquad++) {
/*build the array of the points*/ dirt_grid = DMA_build_DEM(dted,&numvert,num_child,nquad);
/*build the vertex and normal index arrays*/ index_array = DMA_build_index_array(dted,numvert,&num_of_points);
/*build the color index array*/ color_index_array = DMA_build_color_array(dirt_grid,index_array, num_of_points,minelev,maxelev,NUMCOLORS);
/*put together the gset*/ map = pfNewGSet(arena); pfGSetPrimType(map,PFGS_TRISTRIPS); pfGSetAttr(map,PFGS_COLOR4,PFGS_PER_VERTEX,color_array,color_index_array); pfGSetAttr(map,PFGS_COORD3,PFGS_PER_VERTEX,dirt_grid,index_array); // pfGSetAttr(map,PFGS_TEXCOORD2, PFGS_PER_VERTEX, texcoords, texlist); pfGSetNumPrims(map,coll_cnt-1); /* coll_cnt-1 = 75 */ // printf("coll_cnt-1 = %d\n", coll_cnt-1); pfGSetGState(map,gstate);
if(NULL == (prim_cnt_arr = (long *)pfMalloc(sizeof(long)*coll_cnt,arena))){ /* printf("Unable to pfMalloc the prim cnt array\n"); */ exit(0); } for(ix=0;ix<coll_cnt;ix++){ prim_cnt_arr[ix] = (2*row_cnt); } pfGSetPrimLengths(map,prim_cnt_arr);
/*pfDebugPrintGSet(map,1,0);*/ pfAddGSet(geode,map);
} /* End loop over 4 quadrants */
/* set up geode intersect testing mask */ // pfNodeTravMask((pfNode *)geode, PFTRAV_ISECT,DIRT_MASK, PFTRAV_SELF|PFTRAV_DESCEND|PFTRAV_IS_CACHE, PF_SET);
pfAddChild(group[num_child],(pfNode *)geode);
} /* end if num_child > 21) */
if(num_child > 1) { /* printf("Adding group %d to group %d\n", num_child, (num_child+2)/4); */ pfAddChild(group[(num_child + 2)/4], (pfNode *) group[num_child]); }
} /* End loop over num_children */
return(group[1]);
}
ushort *DMA_build_color_array(pfVec3 *dirt_grid,ushort *index_array, int num_of_points,float minelev, float maxelev,int numcolor) /*build the table index for the colors*/ { ushort *ca; int ix,jx; float interval;
if(NULL == (ca = (ushort *)pfMalloc(sizeof(ushort)*num_of_points,arena))){ /* printf("Unable to pfMalloc the color index table \n"); */ exit(0); }
interval = (maxelev-minelev)/numcolor; for(ix=0;ix<num_of_points;ix++){ if(dirt_grid[(index_array[ix])][2] == 0.0f){ /*Sea level*/ ca[ix] = 0; } else{ jx = 1; while(dirt_grid[(index_array[ix])][2] > (minelev + jx*interval)){ jx++; } /* if(jx > numcolor){ printf("Color Table Max elev %20.10f Input %20.10f max table %20.10f\n" ,maxelev,dirt_grid[(index_array[ix])][2],minelev+numcolor*interval); jx = numcolor; } */ ca[ix] = jx; /* printf( "%d (%8.1f %8.1f %8.1f) %4d\n",ix, dirt_grid[ix][0],dirt_grid[ix][1],dirt_grid[ix][2],ca[ix]); */ } } return(ca); }
ushort *DMA_build_index_array(DTED *dted,int numvert,int *num_in_index) /*Build the index array*/ { int coll_cnt,row_cnt; ushort *ia; int cnt=0; int ix,jx;
sscanf(dted->Lat_Count,"%4d",&row_cnt); sscanf(dted->Lng_Count,"%4d",&coll_cnt);
if(NULL == (ia = (ushort *)pfMalloc(sizeof(ushort)*(numvert*2),arena))){ /* printf("Unable to pfMalloc the vertex index table \n"); */ exit(0); }
for(ix=0;ix<coll_cnt;ix++){ for(jx = 0;jx<row_cnt;jx++){ ia[cnt] = ix*row_cnt + jx; cnt++; ia[cnt] = (ix+1)*row_cnt + jx;
/* printf("(%4d,%4d) %6d,%6d (%8.1f %8.1f %8.1f)(%8.1f %8.1f %8.1f)\n",ix,jx, ia[cnt-1],ia[cnt], dirt_grid[(ia[cnt-1])][0],dirt_grid[(ia[cnt-1])][1],dirt_grid[(ia[cnt-1])][2], dirt_grid[(ia[cnt])][0],dirt_grid[(ia[cnt])][1],dirt_grid[(ia[cnt])][2]); */ cnt++; } /*printf("********************END OF STRIP *** %10d\n",cnt); */ } /* printf("Number of vertices %10d (%10hd) Counter value %10d\n",numvert, (ushort)numvert,cnt); */ *num_in_index = cnt; return(ia); }
pfVec3 *DMA_build_normal_array(DTED *dted) /*Build the Normals for shading*/ { return(NULL); }
pfVec4 *DMA_build_color_ramp(int numcolors) /*build the color table, note that 0 is always set to blue*/ { pfVec4 *color_table; int ix;
if(NULL == (color_table = (pfVec4 *)pfMalloc(sizeof(pfVec4)*(numcolors +1),arena))){ /* printf("Unable to pfMalloc the color table \n"); */ exit(0); }
/*make the water*/ color_table[0][0] = color_table[0][1] = 0; color_table[0][2] = color_table[0][3] = 1.0f;
/*set up the rest of the grey scale*/ for (ix=1;ix<=numcolors;ix++){ /* color_table[ix][0] = color_table[ix][1] = color_table[ix][2] = ix * (1.0f/numcolors); */
/* greenish */ /* color_table[ix][0] = 1.1f * ix * (1.0f/numcolors); color_table[ix][1] = 1.5f * ix * (1.0f/numcolors); color_table[ix][2] = 0.9f * ix * (1.0f/numcolors); */ /* brownish */ color_table[ix][0] = 1.5f * ix * (1.0f/numcolors); color_table[ix][1] = 0.9f * ix * (1.0f/numcolors); color_table[ix][2] = 0.6f * ix * (1.0f/numcolors);
color_table[ix][3] = 1.0f; /*printf("Color %3d %10.6f %10.6f %10.6f %10.6f\n",ix,color_table[ix][0],color_table[ix][1],color_table[ix][2],color_table[ix][3]); */ }
return(color_table); }
void DMA_read_header (DTED *dted,int *latitude,int *longitude, int *row_cnt, int *col_cnt, float *lat_spacing,float *long_spacing) /*read the header and compute all of the values*/ { int lng_interval,lat_interval;
/*----- Determine latitude in degrees of DTED grid origin from DSI header */ sscanf(dted->Lat_Origin, "%2d", latitude); sscanf(dted->Lng_Origin, "%3d", longitude); /*----- Determine longitude & latitude spacing of terrain data points */ /*----- from DSI header information. */ /*----- Spacing is given in tenths of arc seconds. */ sscanf(dted->Lng_Interval, "%4d", &lng_interval); sscanf(dted->Lat_Interval, "%4d", &lat_interval);
/*----- Convert terrain data spacing to meters, corrected for latitude */ /*----- assuming a spherical earth model. Longitudinal interval (I) at */ /*----- latitude (L) fractional of equatorial interval (E) as given by */ /*----- I = cos(L) * E . Latitudinal interval is not affected. */ /*----- delta = tenths / (10.0 tenths/s) * (1852.0 m/nm) / (60.0 s/nm) */ *long_spacing = (float)lng_interval / 10.0 * 1852.0 / 60.0 * fcos((*latitude*0.017453293)); *lat_spacing = (float)lat_interval / 10.0 * 1852.0 / 60.0;
sscanf(dted->Lat_Count,"%4d",row_cnt); sscanf(dted->Lng_Count,"%4d",col_cnt);
}
void min_max_elev(DTED *dted, float *minelev, float *maxelev) { int num_col, num_rows; float delta_lng,delta_lat; int ix,jx; int latitude,longitude;
*maxelev = 0.0f; *minelev = 99999999999.9f;
DMA_read_header(dted,&latitude,&longitude,&num_rows,&num_col, &delta_lat,&delta_lng);
for(ix=0;ix <num_col;ix++){ for(jx=0;jx<num_rows;jx++){
if ( (float)dted->terrain[ix][jx] > *maxelev) *maxelev = (float)dted->terrain[ix][jx]; if ( (float)dted->terrain[ix][jx] < *minelev) *minelev = (float)dted->terrain[ix][jx]; } }
}
pfVec3 * DMA_build_DEM(DTED *dted, int *numvert, int num_child, int nquad) /*build the point array from the DTED elevation array*/ {
int num_col, num_row; float delta_lng,delta_lat; pfVec3 *vert_list; int ix,jx, ix_start, ix_end, jx_start, jx_end; int latitude,longitude; int cnt = 0;
DMA_read_header(dted,&latitude,&longitude,&num_row,&num_col, &delta_lat,&delta_lng);
*numvert = num_col*num_row;
/* printf(" %d %d %d\n",num_col,num_row,*numvert); */ if(NULL == (vert_list = (pfVec3 *)pfMalloc(sizeof(pfVec3)*(*numvert),arena))){ /* printf("Unable to pfMalloc the vertex array count %d\n",*numvert); */ exit(0); }
/* for(ix=0;ix <num_col;ix++){ */ /* for(ix=num_child*num_col-num_child;ix <(num_child+1)*num_col-num_child;ix++){ for(jx=0;jx<num_rows;jx++){ */ ix_start = origin[num_child][0]+quadr[nquad][0]*(num_col-1); ix_end = ix_start+num_col; jx_start = origin[num_child][1]+quadr[nquad][1]*(num_row-1); jx_end = jx_start+num_row; /* printf("nquad = %d: ix_start = %d ix_end = %d jx_start = %d jx_end = %d\n", nquad, ix_start,ix_end, jx_start,jx_end); */
for(ix=ix_start;ix < ix_end;ix++){ for(jx=jx_start;jx< jx_end;jx++){
if (dted->terrain[ix][jx] < 0){ /* printf("Value less then zero %d (X %d Z %d)\n", dted->terrain[ix][jx] ,ix,jx); */ dted->terrain[ix][jx] = 0; } vert_list[cnt][0] = ix*delta_lng; vert_list[cnt][1] = jx*delta_lat; vert_list[cnt][2] = (float)dted->terrain[ix][jx]; cnt++; } }
/* printf("Number vertices %d\n",cnt); */ *numvert = cnt;
return(vert_list); }
/**************************************************************************** Filename: dma_dted.c Purpose: Processes DMA DTED Level 1/2 Terrain Data Revised: 90/03/22 Author: Richard M. Prevatt David M. King ........................................................................ Functions: read_DTED_terrain() Reads DTED data using std header read_DTED_normals() Reads vertex normals for DTED data write_DTED_normals() Writes normals for DTED data to file prepare_column_normals()Prepare structures for column normals compute_column_normals()Computes vertex normals fm DTED data DMA_compute_vertex_normal() Computes avg vector of 6 cross prods DMA_crossProduct() Computes cross product of 2 vectors DMA_load_vector_x_z() Load x & z components of vectors DMA_load_vector_y() Load y components of vectors ........................................................................ Usage: include the "dma_dted.h" header file in any files which will reference DTED structures or functions. Link a compiled version of "dma_dted.c" to the applications programs. ........................................................................ Notes: 90/01/15 Created. 90/01/21 Add DTED structure to contain header info & data. Correct skip beyond ACC record. 90/01/22 Add terrain data array dynamic allocation. Complete documentation. 90/01/25 Add functions dealing with terrain normals and cross products. 90/01/28 Use float for all cross products and normals. Complete DMA_compute_vertex_normal() ****************************************************************************/ #include <stdlib.h> #include <stdio.h> #include <math.h> #include "dma_dted.h" /* Constants and Data Structures */
extern void *arena;
/* ======================================================================== */ DTED *DMA_read_DTED_terrain (char dted_path[256], DTED *dted) // char dted_path[256]; /* Fully specified path for datafile*/ // DTED *dted; /* DTED data structure ptr or NULL */ /* ------------------------------------------------------------------------ */ /* Read DMA Digital Terrain Elevation Data */ /* Arguments: fully specified path & filename for the DTED data file. */ /* optional pointer to DTED data structure. NULL causes */ /* this function to dynamically allocate the data structure. */ /* Calls: none. */ /* Returns: pointer to the filled DTED data structure. */ /* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . */ /* Any standard DTED Level 1/2 file can be processed based on information */ /* within the header of the data file. Three files are normally */ /* associated with DTED data: */ /* datafile.head - contains VOL and HDR information. Ignored here. */ /* datafile - contains DSI, ACC and actual terrain data. Used. */ /* datafile.tail - contains VOL and EOF information. Ignored here. */ /* The main data file is processed by this function. Both DSI header info */ /* and terrain elevation data is retained in the DTED data structure. */ /* Accuracy information (ACC) is ignored. */ /* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . */ /* Specification for the DTED data is detailed in: */ /* Defense Mapping Agency Product Specifications for */ /* Digital Terrain Elevation Data (DTED), 2nd Edition, April 1986 */ /* ------------------------------------------------------------------------ */ { FILE *dted_file; /* pointer to data file */ int col, max_col, /* column index and limit variables */ max_row; /* row index and limit variables */
if ( !dted ) { /* No DTED data structure passed in */ dted = (DTED*)pfMalloc(sizeof(DTED),arena); /* Create one */ dted->terrain = NULL; /* Initialize terrain array pointer */ }
/*------------------------------------/* Open DTED data file */
if ( (dted_file = fopen(dted_path,READ)) == NULL ) { fprintf(stderr, "Cannot open DMA DTED data file '%s'\n", dted_path); exit(0); /* Unable... Error and Exit */ } printf("\nReading in %s",dted_path); fflush(stdout);
/*------------------------------------/* Read DTED DSI Record information */ fread(dted, BYTE, 648, dted_file);
fseek(dted_file, 2700, SEEK_CUR); /* Skip Accuracy Description Record */
/*------------------------------------/* Compute max col & row fm DSI info */ sscanf(dted->Lng_Count, "%4d", &max_col); sscanf(dted->Lat_Count, "%4d", &max_row);
/*------------------------------------/* Dynamically allocate terrain array*/ if ( !dted->terrain ) { /* if not already allocated. */ /*----------------------------------/* allocate the column pointers */ dted->terrain = (short**) pfMalloc (max_col * sizeof(short *),arena); for ( col=0; col < max_col; col++ ) /* for each col, allocate storage */ dted->terrain[col] = (short*) pfMalloc (max_row * sizeof(short),arena); }
/*------------------------------------/* For each longitude column */ for ( col=0; col < max_col; col++) { # ifdef TRACE fprintf(stderr,"read_DTED_terrain(): column [%d]\n", col); # endif fseek(dted_file, 8, SEEK_CUR); /* Skip block count fields */ /*----------------------------------/* Read column of data */ fread(dted->terrain[col], BYTE, max_row*2, dted_file); fseek(dted_file, 4, SEEK_CUR); /* Skip checksum field */ }
printf("\nDMA_read_DTED_terrain() # col = %d, # row = %d\n\n", max_col, max_row); fclose(dted_file); /* Close DTED data file */
return(dted); /* Return pointer to DTED structure */
} /* read_DTED_terrain() */
/* ------------------------------------------------------------------------ */
/* ======================================================================== */ DTED *DMA_read_DTED_normals (char norm_path[256],DTED *dted) // char norm_path[256]; /* Fully specified path for datafile*/ // DTED *dted; /* DTED data structure ptr or NULL */ /* ------------------------------------------------------------------------ */ /* Read vertex normals corresponding to DTED terrain data */ /* Arguments: fully specified path & filename for the terrain normals */ /* data file, usually 'datafile.norm' . */ /* optional pointer to DTED data structure. NULL causes */ /* this function to dynamically allocate the data structure. */ /* Calls: none. */ /* Returns: pointer to the filled DTED data structure. */ /* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . */ /* ------------------------------------------------------------------------ */ { FILE *norm_file; /* pointer to data file */ int axis, /* normal component axis {x, y, z} */ col, max_col, /* column index and limit variables */ max_row; /* row index and limit variables */
if ( !dted ) { /* No DTED data structure passed in */ dted = (DTED*)pfMalloc(sizeof(DTED),arena); /* Create one */ dted->normal = NULL; /* Initialize normals array pointer */ }
/*------------------------------------/* Open Terrain Normals data file */ if ( (norm_file = fopen(norm_path,READ)) == NULL ) { fprintf(stderr, "Cannot open DTED Terrain Normals data file '%s'\n", norm_path); exit(-1); /* Unable... Error and Exit */ }
/*------------------------------------/* Read DTED DSI Record information */ fread(dted, BYTE, 648, norm_file);
/*------------------------------------/* Compute max col & row fm DSI info */ sscanf(dted->Lng_Count, "%4d", &max_col); sscanf(dted->Lat_Count, "%4d", &max_row);
/*------------------------------------/* Dynamically allocate normals array*/ if ( !dted->normal ) { /* if not already allocated. */ /*----------------------------------/* allocate the axis pointers */ dted->normal = (float***) pfMalloc ( 3 * sizeof(float**),arena); /*----------------------------------/* for each axis, allocate col ptr */ for ( axis=0; axis < 3; axis++ ) { dted->normal[axis] = (float**) pfMalloc ( max_col * sizeof(float*),arena); /*--------------------------------/* for each col, allocate storage */ for ( col=0; col < max_col; col++ ) dted->normal[axis][col] = (float*) pfMalloc (max_row * sizeof(float),arena); } }
/*------------------------------------/* For each axis {x, y, z} */ for ( axis=0; axis < 3; axis++ ) { /*----------------------------------/* and each longitude column */ for ( col=0; col < max_col; col++) { # ifdef TRACE fprintf(stderr,"read_DTED_normals(): column [%d]\n", col); # endif /*--------------------------------/* Read column of normal components */ fread(dted->normal[axis][col], REAL, max_row, norm_file); } }
fclose(norm_file); /* Close DTED Normals data file */
return(dted); /* Return pointer to DTED structure */
} /* read_DTED_normals() */ /* ------------------------------------------------------------------------ */
/* ======================================================================== */ boolean DMA_write_DTED_normals (char norm_path[256],DTED *dted) // char norm_path[256]; /* Fully specified path for datafile*/ // DTED *dted; /* DTED data structure ptr or NULL */ /* ------------------------------------------------------------------------ */ /* Write vertex normals corresponding to DTED terrain data to file. */ /* Arguments: fully specified path & filename for the terrain normals */ /* data file, usually 'datafile.norm' . */ /* required pointer to DTED data structure. */ /* Calls: none. */ /* Returns: SUCCESS (1) if successful. FAILURE (0) if unsuccessful. */ /* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . */ /* ------------------------------------------------------------------------ */ { FILE *norm_file; /* pointer to data file */ int axis, /* normal component axis {x, y, z} */ col, max_col, /* column index and limit variables */ max_row; /* row index and limit variables */
if ( !dted ) { /* No DTED data structure passed in */ fprintf(stderr, "DTED data structure not specified. \n"); return(FAILURE); /* Error and Return */ }
if ( !dted->normal ) { /* No normals array present */ fprintf(stderr, "DTED normals array not present in data structure. \n"); return(FAILURE); /* Error and Return */ }
/*------------------------------------/* Open Terrain Normals data file */ if ( (norm_file = fopen(norm_path,WRITE)) == NULL ) { fprintf(stderr, "Cannot open DTED Terrain Normals data file '%s'\n", norm_path); return(FAILURE); /* Unable... Error and Return */ }
/*------------------------------------/* Write DTED DSI Record information */ fwrite(dted, BYTE, 648, norm_file);
/*------------------------------------/* Compute max col & row fm DSI info */ sscanf(dted->Lng_Count, "%4d", &max_col); sscanf(dted->Lat_Count, "%4d", &max_row);
/*------------------------------------/* For each axis {x, y, z} */ for ( axis=0; axis < 3; axis++ ) { /*----------------------------------/* and each longitude column */ for ( col=0; col < max_col; col++) { # ifdef TRACE fprintf(stderr,"write_DTED_normals(): axis[%d] column[%d]\n",axis, col); # endif /*--------------------------------/* Write column of normal components */ fwrite(dted->normal[axis][col], REAL, max_row, norm_file); } }
fclose(norm_file); /* Close DTED Normals data file */
return(TRUE); /* Return success of operation */
} /* write_DTED_normals() */ /* ------------------------------------------------------------------------ */
/* ======================================================================== */ boolean DMA_prepare_column_normals (DTED *dted) // DTED * dted;
/* ------------------------------------------------------------------------ */ /* Prepare structures used to calculate column normals corresponding to */ /* DTED terrain data */ /* Arguments: None. */ /* Calls: None. */ /* Returns: SUCCESS (1) if successful. FAILURE (0) if unsuccessful. */ /* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . */ /* ------------------------------------------------------------------------ */ { int axis, /* normal component axis {x, y, z} */ col; /* column index and limit variables */ int max_col,max_row;
if ( !dted ) { /* No DTED data structure passed in */ fprintf(stderr, "DTED data structure not specified. \n"); return(FAILURE); /* Error and Return */ }
if ( !dted->terrain ) { /* No terrain data array present */ fprintf(stderr, "DTED normals array not present in data structure. \n"); return(FAILURE); /* Error and Return */ }
/*------------------------------------/* Compute max col & row fm DSI info */ sscanf(dted->Lng_Count, "%4d", &max_col); sscanf(dted->Lat_Count, "%4d", &max_row);
/*------------------------------------/* Dynamically allocate normals array*/ if ( !dted->normal ) { /* if not already allocated. */ /*----------------------------------/* allocate the axis pointers */ dted->normal = (float***) pfMalloc ( 3 * sizeof(float**),arena); /*----------------------------------/* for each axis, allocate col ptr */ for ( axis=0; axis < 3; axis++ ) { dted->normal[axis] = (float**) pfMalloc ( max_col * sizeof(float*),arena); /*--------------------------------/* for each col, allocate storage */ for ( col=0; col < max_col; col++ ) dted->normal[axis][col] = (float*) pfMalloc (max_row * sizeof(float),arena); } }
return(SUCCESS); /* Return successful completion */
} /* prepare_column_normals() */ /* ------------------------------------------------------------------------ */
/* ======================================================================== */ boolean DMA_compute_column_normals (DTED *dted) // DTED *dted; /* ------------------------------------------------------------------------ */ /* Compute vertex normals for each column using parallel processors. */ /* Arguments: None. */ /* Calls: None. /* Returns: SUCCESS (1) if successful. FAILURE (0) if unsuccessful. */ /* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . */ /* ------------------------------------------------------------------------ */ { int col, row, /* current column and row */ /* column groupings */ /* first and last column of group */ latitude, /* latitude of terrain grid origin */ lng_interval, /* longitudinal data point spacing */ lat_interval; /* latitudinal data point spacing */ int max_row,max_col; boolean v_valid[6]; /* vertex vectors validity array */ double delta_lat, /* terrain data spacing in meters */ delta_lng, v_normal[3], /* resulting vertex normal */ v[6][3]; /* vectors surrounding the vertex */
/*----- Determine latitude in degrees of DTED grid origin from DSI header */
sscanf(dted->Lat_Origin, "%2d", &latitude);
/*----- Determine longitude & latitude spacing of terrain data points */ /*----- from DSI header information. */ /*----- Spacing is given in tenths of arc seconds. */
sscanf(dted->Lng_Interval, "%4d", &lng_interval); sscanf(dted->Lat_Interval, "%4d", &lat_interval);
/*----- Convert terrain data spacing to meters, corrected for latitude */ /*----- assuming a spherical earth model. Longitudinal interval (I) at */ /*----- latitude (L) fractional of equatorial interval (E) as given by */ /*----- I = cos(L) * E . Latitudinal interval is not affected. */ /*----- delta = tenths / (10.0 tenths/s) * (1852.0 m/nm) / (60.0 s/nm) */
delta_lng = (double)lng_interval / 10.0 * 1852.0 / 60.0 * cos(latitude); delta_lat = (double)lat_interval / 10.0 * 1852.0 / 60.0; sscanf(dted->Lat_Count,"%4d",&max_row); sscanf(dted->Lng_Count,"%4d",&max_col);
/*----- Load the x & z components of the surrounding vectors array */
DMA_load_vector_x_z(delta_lat, delta_lng, v);
for ( col=0; col < max_col; col++ ) { /* For each col in group */ for ( row=0; row < max_row; row++ ) {/* and each row in column */
/*----- Load the y component of the surrounding vectors array */
DMA_load_vector_y(dted->terrain, col, row, max_col-1, max_row-1, v, v_valid);
DMA_compute_vertex_normal(v, v_valid, v_normal);
/*----- Save the vertex normal in the dted structure */
dted->normal[0][col][row] = (float) v_normal[0]; dted->normal[1][col][row] = (float) v_normal[1]; dted->normal[2][col][row] = (float) v_normal[2]; } } return(SUCCESS);
} /* compute_column_normals() */ /* ------------------------------------------------------------------------ */
/* ======================================================================== */ int DMA_compute_vertex_normal (double v[6][3],boolean v_flag[6], double result[3] ) // double v[6][3]; /* Input vectors */ // boolean v_flag[6]; /* Validity of each input vector */ // double result[3]; /* Result of average v[i] X ... */ /* ------------------------------------------------------------------------ */ /* Calculate the cross product of up to six vectors, v[1] X v[2] X ... */ /* Arguments: six vectors of three dimensions */ /* six flags showing which vectors contain valid data. */ /* the vector for the resulting cross products averaged. */ /* Calls: DMA_crossProduct() */ /* Returns: SUCCESS (1) if successful. FAILURE (0) if unsuccessful. */ /* ------------------------------------------------------------------------ */ { int ndx1, ndx2; /* Index for vector array */ double cross_product[3], /* Temporary cross_product */ magnitude; /* Magnitude of result vector */
result[0] = 0.0; /* Initialize result to zero */ result[1] = 0.0; result[2] = 0.0;
for ( ndx1=0; ndx1 < 5; ndx1++ ) { /* For each of first 5 vectors */ ndx2 = ndx1 + 1; /* Compute adjacent index */ /*----------------------------------/* if two adjacent vectors are valid */ if ( v_flag[ndx1] && v_flag[ndx2] ) { /*--------------------------------/* Compute the cross product */ DMA_crossProduct( v[ndx1], v[ndx2], cross_product );
result[0] += cross_product[0]; /* Add new cross product to result */ result[1] += cross_product[1]; result[2] += cross_product[2]; } }
/*------------------------------------/* For sixth vector, */ if ( v_flag[5] && v_flag[0] ) { /* if adjacent vectors are valid */ /*----------------------------------/* Compute the cross product */ DMA_crossProduct( v[5], v[0], cross_product );
result[0] += cross_product[0]; /* Add new cross product to result */ result[1] += cross_product[1]; result[2] += cross_product[2]; }
/*------------------------------------/* Compute magnitude */ magnitude = fsqrt( result[0] * result[0] + result[1] * result[1] + result[2] * result[2] );
/*------------------------------------/* Normalize the result vector */ if ( magnitude ) { /* if the magnitude is not zero */ result[0] /= magnitude; result[1] /= magnitude; result[2] /= magnitude; }
return(SUCCESS);
} /* DMA_compute_vertex_normal() */ /* ------------------------------------------------------------------------ */
/* ======================================================================== */ void DMA_crossProduct ( double v1[3],double v2[3],double result[3] ) // double v1[3], v2[3], /* Input vectors */ // result[3]; /* Result of v1 X v2 */ /* ------------------------------------------------------------------------ */ /* Calculate the cross product of two vectors, v1 X v2. */ /* Arguments: two vectors of three dimensions. */ /* the vector for the resulting cross product. */ /* Calls: none. */ /* Returns: none. */ /* ------------------------------------------------------------------------ */ { double magnitude; /* Magnitude of result */
/*------------------------------------/* Compute cross product */ result[0] = v1[1] * v2[2] - v1[2] * v2[1]; result[1] = v1[2] * v2[0] - v1[0] * v2[2]; result[2] = v1[0] * v2[1] - v1[1] * v2[0];
/*------------------------------------/* Compute magnitude */ magnitude = fsqrt( result[0] * result[0] + result[1] * result[1] + result[2] * result[2] );
/*------------------------------------/* Normalize the result vector */ if ( magnitude ) { /* if the magnitude is not zero */ result[0] /= magnitude; result[1] /= magnitude; result[2] /= magnitude; }
} /* DMA_crossProduct() */ /* ------------------------------------------------------------------------ */
/* ======================================================================== */ void DMA_load_vector_x_z ( double delta_x,double delta_z,double v[6][3] ) // double delta_x, /* increment in x direction */ // delta_z, /* increment in z direction */ // v[6][3]; /* vector array */ /* ------------------------------------------------------------------------ */ /* Load array v with x & z components of six vectors based on orientation */ /* and increments in respective directions. */ /* Arguments: increment in x and z direction. */ /* array 'v' to be loaded with surrounding vectors */ /* Calls: none. */ /* Returns: none. */ /* ------------------------------------------------------------------------ */ { v[0][0] = delta_x; /* [2] [1] */ v[0][2] = 0.0; /* * * */ v[1][0] = delta_x; /* | / */ v[1][2] = -delta_z; /* | / */ v[2][0] = delta_x; /* | / */ v[2][2] = -delta_z; /* [3] * - - - * - - - * [0] */ v[3][0] = -delta_x; /* / | */ v[3][2] = 0.0; /* / | */ v[4][0] = -delta_x; /* / | -z */ v[4][2] = delta_z; /* * * ^ */ v[5][0] = 0.0; /* [4] [5] | */ v[5][2] = delta_z; /* +--> x */ } /* ------------------------------------------------------------------------ */
/* ======================================================================== */ void DMA_load_vector_y(short **terrain,int col,int row,int col_max,int row_max, double v[6][3], boolean v_valid[6] ) // short **terrain; /* terrain elev data array */ // int col, row, /* vertex col & row position */ // col_max, row_max; /* maximum col & row in array */ // double v[6][3]; /* vector array */ // boolean v_valid[6]; /* flags to indicate valid vectors */ /* ------------------------------------------------------------------------ */ /* Load array v with the six vectors surrounding the given vertex. */ /* Arguments: pointer to 'terrain' elevation data array */ /* column and row of vertex to be considered */ /* column and row maximum values */ /* array 'v' to be loaded with surrounding vectors */ /* array 'v_valid' to indicate which vectors are valid */ /* Calls: none. */ /* Returns: none. */ /* ------------------------------------------------------------------------ */ { /*------------------------------------/* Determine location of vertex */ if ( col < col_max ) { /* Not in rightmost column */
v[0][1] = (double)(terrain[col+1][row+0] - terrain[col][row]); v_valid[0] = TRUE;
if ( row < row_max ) { /* Not in topmost row */
v[1][1] = (double)(terrain[col+1][row+1] - terrain[col][row]); v[2][1] = (double)(terrain[col+0][row+1] - terrain[col][row]); v_valid[1] = TRUE; v_valid[2] = TRUE;
} else { /* In topmost row */ v_valid[1] = FALSE; v_valid[2] = FALSE; }
} else { /* In rightmost column */
v_valid[0] = FALSE; v_valid[1] = FALSE;
if ( row < row_max ) { /* Not in topmost row */
v[2][1] = (double)(terrain[col+0][row+1] - terrain[col][row]); v_valid[2] = TRUE; }
else /* In topmost row */
v_valid[2] = FALSE; }
if ( col > 0 ) { /* Not in leftmost column */
v[3][1] = (double)(terrain[col-1][row+0] - terrain[col][row]); v_valid[3] = TRUE;
if ( row > 0 ) { /* Not in bottom row */
v[4][1] = (double)(terrain[col-1][row-1] - terrain[col][row]); v[5][1] = (double)(terrain[col+0][row-1] - terrain[col][row]); v_valid[4] = TRUE; v_valid[5] = TRUE;
} else { /* In bottom row */
v_valid[4] = FALSE; v_valid[5] = FALSE; }
} else { /* In leftmost row */
v_valid[3] = FALSE; v_valid[4] = FALSE;
if ( row > 0 ) { /* Not in bottom row */
v[5][1] = (double)(terrain[col+0][row-1] - terrain[col][row]); v_valid[5] = TRUE;
} else /* In bottom row */
v_valid[5] = FALSE; } } /* load_vectors() */ /* ------------------------------------------------------------------------ */
/**************************************************************************** End of File: dma_dted.c lines: 567 ****************************************************************************/
#ifndef _DMA_DTED_H_ #define _DMA_DTED_H_ /**************************************************************************** Filename: dma_dted.h2 Purpose: Processes DMA DTED Level 1/2 Terrain Data Revised: 90/08/09 Author: Richard M. Prevatt ........................................................................ Constants and Structure Declarations for DTED functions ........................................................................ Notes: 90/01/15 Created. 90/01/21 Add DTED structure to contain header info & data. 90/01/22 Add terrain data array dynamic allocation. Complete documentation. 90/01/25 Add pointer for terrain normals array. ****************************************************************************/
/* Constant Declarations ================================================== */ /* #include <stdlib.h> #include <sys/types.h> #include <sys/times.h> #include <ulocks.h> #include <sys/param.h> #include <task.h> #include <sys/prctl.h> #include <pr.h> /* added so that we could use Performer's pfMalloc */
/* Constant Declarations ================================================== */
#define boolean int #define BYTE sizeof(char) #define REAL sizeof(float)
#define SEEK_SET 0 /* Seek from beginning of file */ #define SEEK_CUR 1 /* Seek from current position */ #define SEEK_END 2 /* Seek to end of file */
#define READ "r" #define WRITE "w" #define APPEND "a"
#ifndef FALSE #define FALSE 0 #endif #ifndef TRUE #define TRUE !FALSE #endif #define FAILURE FALSE #define SUCCESS !FAILURE
/* Type Definitions ======================================================= */
/* ======================================================================== */ typedef struct DTED_type { /* DTED data structure type def */ /* ------------------------------------------------------------------------ */ /* DTED is a defined type to provide character storage for DTED Data Set */ /* Identification Record information plus a pointer 'terrain' for the */ /* two dimensional array to contain actual elevation data and a pointer */ /* 'normals' for the three dimensional array to contain terrain normals. */ /* All Data Set Identification (DSI) information is ascii character data. */ /* Terrain elevation data are two byte signed binary numbers (short). */ /* Access to terrain data points has the form : */ /* dted.terrain[col][row] or dted_ptr->terrain[col][row] */ /* Terrain normals are 32 bit floating point numbers (double). */ /* Access to terrain normals has the form : */ /* dted.normal[axis][col][row] or dted_ptr->normal[axis][col][row] */ /* where axis = 0, 1, 2 for x, y, z respectively. */ /* . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . */ /* Specification for the DTED data is detailed in: */ /* Defense Mapping Agency Product Specifications for */ /* Digital Terrain Elevation Data (DTED), 2nd Edition, April 1986 */ /* ------------------------------------------------------------------------ */ char Sentinel[3], /* 'DSI' Recognition Sentinel */ Class[1], /* Security Classification Code */ SecMark[2], /* Security Control & Release Mark */ SecDesc[27], /* Security Handling Description */ rsvd_1[26], /* Reserved for future use */ Level[5], /* DMA Series Designator for level */ RefNo[15], /* Unique reference number */ rsvd_2[8], /* Reserved for future use */ Edition[2], /* Data Edition Number */ MergeVersion[1], /* Match/Merge Version */ MaintDate[4], /* Maintenance Date {YYMM} */ MergeDate[4], /* Match/Merge Date {YYMM} */ MaintDesc[4], /* Maintenance Description Code */ Producer[8], /* Producer Code */ rsvd_3[16], /* Reserved for future use */ ProdSpecNo[9], /* Product Specification Stock No */ ProdSpecChg[2], /* Product Specification Change No */ ProdSpecDate[4], /* Product Specification Date {YYMM} */ VertDatum[3], /* Vertical Datum {MSL most cases} */ HorizDatum[5], /* Horizontal Datum Code */ CollectSys[10], /* Digitizing Collection System */ CompileDate[4], /* Compilation Date {YYMM} */ rsvd_4[22], /* Reserved for future use */ Lat_Origin[9], /* Latitude of origin of data */ Lng_Origin[10], /* Longitude of origin of data */ Lat_SW[7], /* Latitude of SW corner of data */ Lng_SW[8], /* Longitude of SW corner of data */ Lat_NW[7], /* Latitude of NW corner of data */ Lng_NW[8], /* Longitude of NW corner of data */ Lat_NE[7], /* Latitude of NE corner of data */ Lng_NE[8], /* Longitude of NE corner of data */ Lat_SE[7], /* Latitude of SE corner of data */ Lng_SE[8], /* Longitude of SE corner of data */ Orientation[9], /* Clockwise orientation angle fm N */ Lat_Interval[4], /* Latitude interval {tenth secs} */ Lng_Interval[4], /* Longitude interval {tenth secs} */ Lat_Count[4], /* Number of latitude lines {rows} */ Lng_Count[4], /* Number of longitude lines {cols} */ CellCover[2], /* Partial cell indicator */ rsvd_5[357]; /* Reserved for future use */ short **terrain; /* Pointer for Terrain Elev Data */ float ***normal; /* Pointer for Terrain normals */ } DTED; /* defined type name */ /* ------------------------------------------------------------------------ */
/* Global Shared Data Structures ========================================== */
#ifdef _SHOULD_NOT_ extern DTED *dted; /* DTED data structure ptr or NULL */ extern int max_col, max_row; /* column and row limit variables */ #endif /* ------------------------------------------------------------------------ */
/* Function Declarations ================================================== */
#ifdef __cplusplus extern "C" { #endif
DTED *DMA_read_DTED_terrain ( char[], DTED* ); DTED *DMA_read_DTED_normals ( char[], DTED* ); boolean DMA_write_DTED_normals ( char[], DTED* ); boolean DMA_prepare_column_normals ( ); boolean DMA_compute_column_normals ( ); int DMA_compute_vertex_normal ( double[6][3], boolean[6], double[3] ); void DMA_load_vector_x_z ( double, double, double[6][3] ); void DMA_load_vector_y ( short**, int, int, int, int, double[6][3], boolean[6] ); void DMA_crossProduct ( double[3], double[3], double[3] );
#ifdef __cplusplus } #endif
#ifdef _DMA_LOADER_ /* Loader functions */ pfVec4 *DMA_build_color_ramp(int); void DMA_read_header (DTED *,int *,int *, int *,int *, float *, float *); void min_max_elev (DTED *, float *, float *); pfVec3 *DMA_build_DEM(DTED *, int *,int, int); pfVec3 *DMA_build_normal_array(DTED *); ushort *DMA_build_color_array(pfVec3 *,ushort *,int ,float ,float,int); ushort *DMA_build_index_array(DTED *,int,int*);
#endif
/**************************************************************************** End of File: dma_dted.h lines: 117 ****************************************************************************/
#endif
extern "C" { void *pfMalloc(size_t nbytes,void *arena); }
This archive was generated by hypermail 2.0b2 on Mon Aug 10 1998 - 17:51:07 PDT