Writing LAS file in own C++ program

1,802 views
Skip to first unread message

George Vosselman

unread,
Sep 24, 2013, 3:52:12 PM9/24/13
to last...@googlegroups.com
Hi Martin,
 
As discussed in Utrecht I want to incorporate reading/writing of LAS files in our own point processing libraries.
 
I've copied a few lines of code from txt2las.cpp and other laslib files (code below). When I try to open the written laz file with lasview I get the bounding box, but points are not shown. However, when I make a copy of the laz file with
 
las2las -i test.laz -o test2.laz
 
the output of las2las does show the data I've converted.
I've been printing the headers of the file created by my program and the one created by las2las and can't find a difference.
 
Any clue what may go wrong?
 
Cheers,
 
George.
 
 
-----------
 
        #include "LaserPoints.h"  // Our own class definition
        #include "laswriter.hpp"
       
        LASpoint               laspoint;
        LASwriteOpener   laswriteopener;
        LASheader           lasheader;
        LaserPoints          point_set;
        LaserPoints::iterator point;
 
        point_set.Read("test,laser"); // Input of points in own data structure
 
        // Create las file header
        lasheader.clean();
        lasheader.point_data_format = 0; // only X Y Z
        lasheader.point_data_record_length = 20;
        lasheader.number_of_point_records = point_set->size();
        lasheader.number_of_points_by_return[0] = point_set->size(); // All first echo
        lasheader.set_bounding_box(point_set->DataBounds().Minimum().X(),
                             point_set->DataBounds().Minimum().Y(),
                             point_set->DataBounds().Minimum().Z(),
                             point_set->DataBounds().Maximum().X(),
                             point_set->DataBounds().Maximum().Y(),
                             point_set->DataBounds().Maximum().Z());
        
        // Open the laz file
        laswriteopener.set_file_name("test.laz");
        LASwriter* laswriter = laswriteopener.open(&lasheader);
        if (laswriter == 0) {
          fprintf(stderr, "ERROR: could not open laswriter\n");
          exit(0);
        }
        // Write all points
        laspoint.init(&lasheader, lasheader.point_data_format,
                lasheader.point_data_record_length, &lasheader);
                
        for (point=point_set->begin(); point!=point_set->end(); point++) {
          laspoint.x = (int) (100.0 * point->X() + 0.001);
          laspoint.y = (int) (100.0 * point->Y() + 0.001);
          laspoint.z = (int) (100.0 * point->Z() + 0.001);
         
          laswriter->write_point(&laspoint);
        }
        
        // Close laz file
        laswriter->close();
 
 
 

 

 

Martin Isenburg

unread,
Sep 26, 2013, 7:00:36 AM9/26/13
to LAStools - efficient command line tools for LIDAR processing
Hello George,

wonderful little bug. (-: And tricky to find. I attach the lasinfo
output of the small sample that you send me to the end of this email.
When you call

> lasheader.set_bounding_box(point_set->DataBounds().Minimum().X(),
> point_set->DataBounds().Minimum().Y(),
> point_set->DataBounds().Minimum().Z(),
> point_set->DataBounds().Maximum().X(),
> point_set->DataBounds().Maximum().Y(),
> point_set->DataBounds().Maximum().Z());

then internally the offset x y z is set. In the example it is set to 0
400000 0. That means that all Y coordinates will be translated by
400000 implicitly. Until now everything is correct. But then you
manually do the conversion to integers:

> laspoint.x = (int) (100.0 * point->X() + 0.001);
> laspoint.y = (int) (100.0 * point->Y() + 0.001);
> laspoint.z = (int) (100.0 * point->Z() + 0.001);

And now your Y coordinates are off by 400000 because the offset would
need to be subtracted before the multiply. There is automated code for
setting the integer X, Y, and Z of each point. So just replace the
three lines above with

laspoint.set_x(point->X());
laspoint.set_y(point->Y());
laspoint.set_z(point->Z());

Cheers,

Martin

===============================================

D:\lastools\bin>lasinfo own_code_with_laslib.laz
reporting all LAS header entries:
file signature: 'LASF'
file source ID: 0
global_encoding: 0
project ID GUID data 1-4: 00000000-0000-0000-0000-000000000000
version major.minor: 1.2
system identifier: 'LAStools (c) by Martin Isenburg'
generating software: 'laser2laz'
file creation day/year: 0/0
header size: 227
offset to point data: 227
number var. length records: 0
point data format: 2
point data record length: 26
number of point records: 4163
number of points by return: 4163 0 0 0 0
scale factor x y z: 0.01 0.01 0.01
offset x y z: 0 400000 0
min x y z: 93347.21 436270.26 2.85
max x y z: 93358.27 436279.17 19.69
LASzip compression (version 2.2r0 c2 50000): POINT10 2 RGB12 2
reporting minimum and maximum for all LAS point record entries ...
X 9334721 9335827
Y 43627026 43627917
Z 285 1969
intensity 0 0
edge_of_flight_line 0 0
scan_direction_flag 0 0
number_of_returns_of_given_pulse 1 1
return_number 1 1
classification 0 0
scan_angle_rank 0 0
user_data 0 0
point_source_ID 0 0
Color R 17 255
G 33 255
B 26 255
WARNING: 4163 points outside of header bounding box
overview over number of returns of given pulse: 4163 0 0 0 0 0 0
histogram of classification of points:
4163 Created, never classified (0)
real max y larger than header max y by 400000.000000
> --
> Download LAStools at
> http://lastools.org
> http://rapidlasso.com
> Be social with LAStools at
> http://facebook.com/LAStools
> http://twitter.com/LAStools
> http://linkedin.com/groups/LAStools-4408378
> Manage your settings at
> http://groups.google.com/group/lastools/subscribe
Reply all
Reply to author
Forward
0 new messages