txt2las.cpp

Lidar数据处理时

/*
===============================================================================

  FILE:  txt2las.cpp
  
  CONTENTS:
  
		This tool parses LIDAR data as it is typically stored in standard ASCII
		formats and converts it into the more efficient binary LAS format. The
		tool operates in two passes. The first pass counts the points, measures
		their bounding box, and - if applicable - creates the histogram for the
		number of returns.

  PROGRAMMERS:
  
    martin [e-mail]
  
  COPYRIGHT:
  
    copyright (C) 2007  martin [e-mail]
    
    This software is distributed in the hope that it will be useful,
    but WITHOUT ANY WARRANTY; without even the implied warranty of
    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
  
  CHANGE HISTORY:
  
    13 July 2007 -- single pass if output is to file by using fopen("rb+" ... 
    25 June 2007 -- added warning in case that quantization causes a sign flip
    13 June 2007 -- added 'e' and 'd' for the parse string
    26 February 2007 -- created sitting in the SFO lounge waiting for LH 455
  
===============================================================================
*/

#include <stdio.h>
#include <stdlib.h>
#include <string.h>

#include "laswriter.h"

void usage()
{
  fprintf(stderr,"usage:\n");
  fprintf(stderr,"txt2las -parse tsxyz lidar.txt\n");
  fprintf(stderr,"txt2las -parse txyzar lidar.txt.gz lidar.laz\n");
  fprintf(stderr,"txt2las -parse xyz -scale 0.02 -i lidar.txt -o lidar.laz\n");
  fprintf(stderr,"txt2las -parse xyzsst -verbose -xyz_scale 0.05 lidar.txt\n");
  fprintf(stderr,"txt2las -parse xsysz -xyz_scale 0.02 0.02 0.01 lidar.txt\n");
  fprintf(stderr,"txt2las -h\n");
  fprintf(stderr,"---------------------------------------------\n");
  fprintf(stderr,"The '-parse tsxyz' flag specifies how to interpret\n");
  fprintf(stderr,"each line of the ASCII file. For example, 'tsxyzssa'\n");
  fprintf(stderr,"means that the first number is the gpstime, the next\n");
  fprintf(stderr,"number should be skipped, the next three numbers are\n");
	fprintf(stderr,"the x, y, and z coordinate, the next two should be\n");
	fprintf(stderr,"skipped, and the next number is the scan angle.\n");
	fprintf(stderr,"The other supported entries are i - intensity,\n");
	fprintf(stderr,"n - number of returns of given pulse, r - number\n");
	fprintf(stderr,"of return, c - classification, u - user data, and\n");
	fprintf(stderr,"p - point source ID, e - edge of flight line flag, and\n");
	fprintf(stderr,"d - direction of scan flag.\n");
  fprintf(stderr,"---------------------------------------------\n");
  fprintf(stderr,"The '-scale 0.02' flag specifies the quantization. The\n");
  fprintf(stderr,"default value of 0.01 means that the smallest increment\n");
  fprintf(stderr,"two between coordinates is 0.01. If measurements are in\n");
  fprintf(stderr,"meters this corresponds to centimeter accuracy, which is\n");
  fprintf(stderr,"commonly considered sufficient for LIDAR data.\n");
  fprintf(stderr,"---------------------------------------------\n");
  fprintf(stderr,"Other parameters such as '-xyz_offset 500000 2000000 0'\n");
  fprintf(stderr,"or '-xyz_scale 0.02 0.02 0.01' or '-file_creation 67 2003'\n");
  fprintf(stderr,"or '-system_identifier \"Airborne One Leica 50,000 Hz\"'\n");
  fprintf(stderr,"or '-generating_software \"TerraScan\"' are available too.\n");
  fprintf(stderr,"Read the source code for documentation.\n");
  exit(1);
}

static inline void VecUpdateMinMax3dv(double min[3], double max[3], const double v[3])
{
  if (v[0]<min[0]) min[0]=v[0]; else if (v[0]>max[0]) max[0]=v[0];
  if (v[1]<min[1]) min[1]=v[1]; else if (v[1]>max[1]) max[1]=v[1];
  if (v[2]<min[2]) min[2]=v[2]; else if (v[2]>max[2]) max[2]=v[2];
}

static inline void VecCopy3dv(double v[3], const double a[3])
{
  v[0] = a[0];
  v[1] = a[1];
  v[2] = a[2];
}

static bool parse(const char* parse_string, const char* line, double* xyz, LASpoint* point, double* gps_time)
{
	int temp_i;
	float temp_f;
	const char* p = parse_string;
	const char* l = line;

	while (p[0])
	{
		if (p[0] == 'x') // we expect the x coordinate
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%lf", &(xyz[0])) != 1) return false;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'y') // we expect the y coordinate
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%lf", &(xyz[1])) != 1) return false;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'z') // we expect the x coordinate
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%lf", &(xyz[2])) != 1) return false;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 's') // we expect a string or a number that we don't care about
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'i') // we expect the intensity
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%f", &temp_f) != 1) return false;
			if (temp_f < 0.0f || temp_f > 65535.0f) fprintf(stderr, "WARNING: intensity %g is out of range of unsigned short\n", temp_f);
			point->intensity = (unsigned short)temp_f;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'a') // we expect the scan angle
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%f", &temp_f) != 1) return false;
			if (temp_f < -128.0f || temp_f > 127.0f) fprintf(stderr, "WARNING: scan angle %g is out of range of char\n", temp_f);
			point->scan_angle_rank = (char)temp_f;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'n') // we expect the number of returns of given pulse
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%d", &temp_i) != 1) return false;
			if (temp_i < 0 || temp_i > 7) fprintf(stderr, "WARNING: return number %d is out of range of three bits\n", temp_i);
			point->number_of_returns_of_given_pulse = temp_i & 7;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'r') // we expect the number of the return
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%d", &temp_i) != 1) return false;
			if (temp_i < 0 || temp_i > 7) fprintf(stderr, "WARNING: return number %d is out of range of three bits\n", temp_i);
			point->return_number = temp_i & 7;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'c') // we expect the classification
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%d", &temp_i) != 1) return false;
			if (temp_i < 0 || temp_i > 255) fprintf(stderr, "WARNING: classification %d is out of range of unsigned char\n", temp_i);
			point->classification = (unsigned char)temp_i;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'u') // we expect the user data
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%d", &temp_i) != 1) return false;
			if (temp_i < 0 || temp_i > 255) fprintf(stderr, "WARNING: user data %d is out of range of unsigned char\n", temp_i);
			point->user_data = temp_i & 255;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'p') // we expect the point source ID
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%d", &temp_i) != 1) return false;
			if (temp_i < 0 || temp_i > 65535) fprintf(stderr, "WARNING: point source ID %d is out of range of unsigned short\n", temp_i);
			point->point_source_ID = temp_i & 65535;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'e') // we expect the edge of flight line flag
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%d", &temp_i) != 1) return false;
			if (temp_i < 0 || temp_i > 1) fprintf(stderr, "WARNING: edge of flight line flag %d is out of range of boolean flag\n", temp_i);
			point->edge_of_flight_line = (temp_i ? 1 : 0);
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 'd') // we expect the direction of scan flag
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%d", &temp_i) != 1) return false;
			if (temp_i < 0 || temp_i > 1) fprintf(stderr, "WARNING: direction of scan flag %d is out of range of boolean flag\n", temp_i);
			point->scan_direction_flag = (temp_i ? 1 : 0);
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else if (p[0] == 't') // we expect the gps time
		{
			while (l[0] && (l[0] == ' ' || l[0] == ',' || l[0] == '\t')) l++; // first skip white spaces
			if (l[0] == 0) return false;
			if (sscanf(l, "%lf", gps_time) != 1) return false;
			while (l[0] && l[0] != ' ' && l[0] != ',' && l[0] != '\t') l++; // then advance to next white space
		}
		else
		{
	    fprintf(stderr, "ERROR: next symbol '%s' unknown in parse control string\n", p);
		}
		p++;
	}
	return true;
}

#ifdef _WIN32
extern "C" FILE* fopenGzipped(const char* filename, const char* mode);
#endif

int main(int argc, char *argv[])
{
  int i;
  bool dry = false;
  bool verbose = false;
	char* file_name_in = 0;
	char* file_name_out = 0;
	bool itxt = false;
	bool olas = false;
	bool olaz = false;
	double xyz_min[3];
	double xyz_max[3];
	double xyz_scale[3] = {0.01,0.01,0.01};
	double xyz_offset[3] = {0.0,0.0,0.0};
  unsigned int number_of_point_records = 0;
  unsigned int number_of_points_by_return[8] = {0,0,0,0,0,0,0,0};
	char* parse_string = "xyz";
	int file_creation_day = 0;
	int file_creation_year = 0;
	char* system_identifier = 0;
	char* generating_software = 0;
#define MAX_CHARACTERS_PER_LINE 512
  char line[MAX_CHARACTERS_PER_LINE];
	double xyz[3];
	LASpoint point;
	double gps_time;

  for (i = 1; i < argc; i++)
  {
    if (strcmp(argv[i],"-dry") == 0)
    {
      dry = true;
    }
    else if (strcmp(argv[i],"-verbose") == 0)
    {
      verbose = true;
    }
    else if (strcmp(argv[i],"-h") == 0)
    {
      usage();
    }
    else if (strcmp(argv[i],"-parse") == 0)
    {
			i++;
      parse_string = argv[i];
    }
    else if (strcmp(argv[i],"-scale") == 0)
    {
			i++;
			sscanf(argv[i], "%lf", &(xyz_scale[2]));
 			xyz_scale[0] = xyz_scale[1] = xyz_scale[2];
    }
    else if (strcmp(argv[i],"-xyz_scale") == 0)
    {
			i++;
			sscanf(argv[i], "%lf", &(xyz_scale[0]));
			i++;
			sscanf(argv[i], "%lf", &(xyz_scale[1]));
			i++;
			sscanf(argv[i], "%lf", &(xyz_scale[2]));
    }
    else if (strcmp(argv[i],"-xyz_offset") == 0)
    {
			i++;
			sscanf(argv[i], "%lf", &(xyz_offset[0]));
			i++;
			sscanf(argv[i], "%lf", &(xyz_offset[1]));
			i++;
			sscanf(argv[i], "%lf", &(xyz_offset[2]));
    }
    else if (strcmp(argv[i],"-i") == 0)
    {
      i++;
      file_name_in = argv[i];
    }
    else if (strcmp(argv[i],"-itxt") == 0)
    {
      itxt = true;
    }
    else if (strcmp(argv[i],"-o") == 0)
    {
      i++;
      file_name_out = argv[i];
    }
    else if (strcmp(argv[i],"-olas") == 0)
    {
      olas = true;
    }
    else if (strcmp(argv[i],"-olaz") == 0)
    {
      olaz = true;
    }
    else if (strcmp(argv[i],"-file_creation") == 0)
    {
			i++;
			sscanf(argv[i], "%d", &file_creation_day);
			i++;
			sscanf(argv[i], "%d", &file_creation_year);
		}
    else if (strcmp(argv[i],"-system_identifier") == 0 || strcmp(argv[i],"-sys_id") == 0)
    {
			i++;
			system_identifier = argv[i];
		}
    else if (strcmp(argv[i],"-generating_software") == 0 || strcmp(argv[i],"-gen_soft") == 0)
    {
			i++;
			generating_software = argv[i];
		}
		else if (i == argc - 2 && file_name_in == 0 && file_name_out == 0)
		{
			file_name_in = argv[i];
		}
		else if (i == argc - 1 && file_name_in == 0 && file_name_out == 0)
		{
			file_name_in = argv[i];
		}
		else if (i == argc - 1 && file_name_in && file_name_out == 0)
		{
			file_name_out = argv[i];
		}
	}

	// create output file name if none specified and no piped output requested (i.e. neither -olas nor -olaz)

	if (file_name_out == 0 && !olas && !olaz)
	{
		int len = strlen(file_name_in);
		file_name_out = strdup(file_name_in);
		if (file_name_out[len-3] == '.' || file_name_out[len-2] == 'g' || file_name_out[len-1] == 'z')
		{
			len = len - 4;
		}
		while (len > 0 && file_name_out[len] != '.')
		{
			len--;
		}
		file_name_out[len] = '.';
		file_name_out[len+1] = 'l';
		file_name_out[len+2] = 'a';
		file_name_out[len+3] = 's';
		file_name_out[len+4] = '\0';
	}

  // make sure that input and output are not *both* piped

  if (file_name_in == 0 && file_name_out == 0)
  {
		fprintf(stderr, "ERROR: input and output cannot both be pipes\n");
		exit(1);
  }

  // here we make *ONE* big switch. that replicates some code but makes it simple.

  if (file_name_out == 0)
  {
    // because the output goes to a pipe we have to precompute the header
    // information with an additional pass. the input must be a file.

    // open input file for first pass

	  FILE* file_in;

		if (strstr(file_name_in, ".gz"))
		{
#ifdef _WIN32
			file_in = fopenGzipped(file_name_in, "r");
#else
			fprintf(stderr, "ERROR: no support for gzipped input\n");
			exit(1);
#endif
		}
		else
		{
			file_in = fopen(file_name_in, "r");
		}

    if (file_in == 0)
		{
			fprintf(stderr, "ERROR: could not open '%s' for first pass\n",file_name_in);
			exit(1);
		}

	  // create a cheaper parse string that only looks for 'x' 'y' 'z' and 'r'

	  char* parse_less = strdup(parse_string);
	  for (i = 0; i < (int)strlen(parse_string); i++)
	  {
		  if (parse_less[i] != 'x' && parse_less[i] != 'y' && parse_less[i] != 'z' && parse_less[i] != 'r') 
		  {
			  parse_less[i] = 's';
		  }
	  }
	  do
	  {
		  parse_less[i] = '\0';
		  i--;
	  } while (parse_less[i] == 's');

	  // first pass to figure out the bounding box and number of returns

    fprintf(stderr, "first pass over file '%s' with parse '%s'\n", file_name_in, parse_less);

	  // read the first line

	  while (fgets(line, sizeof(char) * MAX_CHARACTERS_PER_LINE, file_in))
	  {
		  if (parse(parse_less, line, xyz, &point, &gps_time))
		  {
			  // init the bounding box
			  VecCopy3dv(xyz_min, xyz);
			  VecCopy3dv(xyz_max, xyz);