Heindl Solutions - (Fast) parsing of a point cloud text file

Point cloud

HALCON offers versatile tools for processing 3D point clouds. In recent months, the number of projects that we implemented and that involved 3D point clouds increased by a large factor. In this post, we would like to tackle one of the first problems one might have to solve when using point clouds in HALCON: How to import text-based point cloud files. This is no problem at all for many common point cloud formats which can be read with the operator read_object_model_3d, for example off, ply, or obj. But in many cases you will received text files in a custom text format. For this article, we will process such a text file manually.

The text file

For a proof-of-concept project, we cooperated with Wenglor to acquire 3D scans of a tire.

The specifications of the Wenglor 3D profile sensor are:

Used profile sensor	`MLWL222`
Working distance	40cm / 16"
Light source	Red Laser
Speed of conveyor belt:	200mm/s
Resolution Z	around 70mm
Resolution X	around 180um

From one of the scans got a 65MB text file containing a point cloud of a tire in the following layout:

-151.7266   -633.6000   -25.2167   0006
-151.8364   -633.6000   -24.8572   0006
-152.0389   -633.6000   -24.8338   0005
-158.3584   -633.0000   -26.4589   0013
-130.1936   -632.7000   -25.0663   0009
-130.4456   -632.7000   -25.2497   0012
-130.6545   -632.7000   -25.2521   0011
-140.1130   -632.7000   -25.5433   0013
-140.3224   -632.7000   -25.5458   0012
⋮

Each line in this text file represents a point. Columns 1, 2, and 3 are the point coordinates X, Y, Z. The 4th column is the associated gray value.

In the following, we will give a basic example on how to read and parse this text file. The first attempt would have taken multiple hours to read the whole file. We will improve the reading and parsing speed in multiple steps and describe the various problems that have to be solved. After all optimizations, we get a reading time of 21s, directly in HDevelop.

A first attempt

read_point_cloud_simple_v0.hdev

* read_point_cloud_simple.hdev, Version 0
Filename := '../tire-heindl-solutions.txt'
count_seconds (Then) 1
open_file (Filename, 'input', FileHandle) 2
X := [] 3
Y := []
Z := []
Gray := []
while (true)
    fread_string (FileHandle, String, IsEOF)
    X := [X, number(String)] 4
    fread_string (FileHandle, String, IsEOF)
    Y := [Y, number(String)]
    fread_string (FileHandle, String, IsEOF)
    Z := [Z, number(String)]
    fread_string (FileHandle, String, IsEOF)
    Gray := [Gray, number(String)]
    if (IsEOF)
        break
    endif
endwhile
close_file (FileHandle)
count_seconds (Now)
DurationSeconds := Now - Then
Msg := 'opening file ' + Filename + ' in ' + DurationSeconds + 's'
dev_inspect_ctrl (Msg) 5

* runtime on machine 'MEGA'
* aborted after 1690s after reading first 75000 entries only (4.4% of 1700000 total)

In this and the following improved scripts, the runtime is measured as follows:
Count number of seconds since a certain fixed point in time (Then) ➡ perform task ➡ count number of seconds since a certain fixed point in time (Now). The duration of the task is then DurationSeconds := Now - Then

A common pattern for reading a file in HALCON is this:

open_file (Filename, 'input', FileHandle)
while (true)
    fread_string (FileHandle, String, IsEOF)
    * parse/store `String`
    *   or
    fread_line (FileHandle, Line, IsEOF)
    * parse/store `String`
    if (IsEOF)
        break
    endif
endwhile
close_file (FileHandle)

The control variables X, Y, Z shall receive collect the parsed X, Y, and Z coordinates of all points. Gray shall collect the parsed gray values between 0 and 255.

Note the conversion from a string to a number. It will become clear that this conversion is flawed.

It is a helpful pattern to store some logging output in a control variable Msg and inspect this variable with dev_inspect_ctrl afterwards.

The runtime of this first attempt to read all the entries would be multiple hours. We aborted this script after about 30 minutes. The first step in optimizing the runtime of the script is to disable updating the program counter, the variable window, and the automatic object output to the active window using dev_update_off

First step: Disable HDevelop updates during script execution

read_point_cloud_simple_v1.hdev

* read_point_cloud_simple.hdev, Version 1
dev_update_off ()

Filename := '../tire-heindl-solutions.txt'
⋮

* runtime on machine 'MEGA'
* interrupted after 200 s, read only first 75000 entries (4.4 % of 1700000 total)

Now the same first 75000 entries can be read in 200s on my test machine. This is much faster, but still way too slow.

Fixing wrong number conversion

After the first few iterations, the result variables X, Y, Z, Gray contain the following entries:

X = [-151.7266, -151.8364, -152.0389, -158.3584, -130.1936, -130.4456]
Y = [-633.6, -633.6, -633.6, -633.0, -632.7, -632.7]
Z = [-633.6, -633.6, -633.6, -633.0, -632.7, -632.7, -25.2497]
Gray = [6, 6, 5, 11, '0009', 10]

Compare this to the input text file: tire-heindl-solutions.txt

-151.7266   -633.6000   -25.2167   0006
-151.8364   -633.6000   -24.8572   0006
-152.0389   -633.6000   -24.8338   0005
-158.3584   -633.0000   -26.4589   0013
-130.1936   -632.7000   -25.0663   0009
-130.4456   -632.7000   -25.2497   0012
⋮

Most entries seem to be OK, but Gray[5] = '0009' is definitely not a number. Even worse, Gray[4] = 11 was parsed as integer number but would be expected to be 13 by most users.

The reason for this is, that the HDevelop operation number parses strings starting with 0 as octal numbers ('001' -> 1, '002' -> 2, '007' -> 7, '008' -> '008', '009' -> '009', '010' -> 8, '011' -> 9, …). This behavior is especially dangerous: For example a number '010', … will be converted to 8 which is often unexpected and hard to detect in later execution steps. number (or the operator tuple_number) should thus be fixed with:

* trim white-space left and right
* trim new line character
* trim leading zeros (which would indicate an octal number to 'tuple_number')
Number := number(regexp_replace(String,'^\\s*0*(.+?)\\s*\n*$', '$1'))

Performance problems

One might expect the overhead of the single fread_string operators to accumulate, as each call to fread_string returns a few bytes per call only. Now we expected that each call to fread_string results in a call to the Windows kernel's ReadFile function.

We tested our assumption with Visual Studio Debugger:

Click to see the details on how to set a break point on kernel32.dll's ReadFile

Launch HDevelop, open read_point_
Execute the script stepwise (F6) until the program counter (green arrow) is at the first call of the operator fread_string
Open Visual Studio
Attach to process: HDevelop
Debug / Break All
Add a new function breakpoint: {,,kernel32.dll}ReadFile
Debug / Continue
In HDevelop: Execute the fread_string operator
Now the breakpoint is hit and Visual Studio shows the disassembly of ReadFile. Having a look at the CPU register BX, the buffer size used by HALCON on x64-win64 seems to be 0x1000 (=4096) bytes per call to ReadFile, thus a call to fread_string results in a call to the Windows kernel function ReadFile approx. every 400th call to fread_string.

Now, our first though was that too many of these calls to ReadFile would be expensive because with each call a trip from user-mode to kernel-mode and back is necessary (ReadFile is a kernel function!). But the HALCON library seems to use an internal buffer or 4096 bytes, which is a good idea. After this tests, we conclude that the fread_string implementation does not seem to cause to much performance overhead. ✅

The more important reason for performance bottlenecks in the script above are statements in HDevelop in the form X := [X, Number[0]] . The internal memory required to store an additional entry in X has to be increased with every call to this statement, and that is an expensive operation. We recommend to pre-allocate an array which is large enough before reading the entries and setting individual tuple entries with X[Index0] := Number[0]. One must then not forget to shrink the probably too large pre-allocated array to the correct number of actually read entries. With this optimization, the runtime to read all entries is 176 seconds:

read_point_cloud_simple_v4.hdev

dev_update_off ()

Filename := '../tire-heindl-solutions.txt'
count_seconds (Then)
open_file (Filename, 'input', FileHandle)

* choose LargeNum at least the size of the expected number of coordinates,
* thus, no rellocations are necessary during insertions into the tuples
* X, Y, Z, Gray. Expected memory allocated: LargeNum * 4 * 8 bytes = 305 MB
LargeNum := 10000000
X := gen_tuple_const(LargeNum, 0.0)
Y := gen_tuple_const(LargeNum, 0.0)
Z := gen_tuple_const(LargeNum, 0.0)
Gray := gen_tuple_const(LargeNum, 0)
Index0 := 0 3
while (true)
    fread_line (FileHandle, OutLine, IsEOF) 1
    if (IsEOF)
        break
    endif
    SplitString := split(OutLine, ' ') 1
    * trim white-space left and right
    * trim new line character
    * trim leading zeros (which would indicate an octal number to 'tuple_number')
    Number := number(regexp_replace(SplitString,'^\\s*0*(.+?)\\s*\n*$', '$1'))
    
    * check if coordinates/gray values are reasonable:
    IndexIsString := find(type_elem(Number),H_TYPE_STRING) 2
    if (IndexIsString > -1)
        throw ('could not convert "' + SplitString[IndexIsString] + '" to a number')
    endif
    
    X[Index0] := Number[0]
    Y[Index0] := Number[1]
    Z[Index0] := Number[2]
    Gray[Index0] := Number[3]
    Index0 := Index0 + 1 3
endwhile
* remember that we initialized the tuples probably too large, we have
* to shrink-to-fit these:
X := X[0:Index0-1]
Y := Y[0:Index0-1]
Z := Z[0:Index0-1]
Gray := Gray[0:Index0-1]
close_file (FileHandle)
count_seconds (Now)
DurationSeconds := Now - Then
Msg := 'opening file ' + Filename + ' in ' + DurationSeconds + 's'
dev_inspect_ctrl (Msg)

* runtime on machine 'MEGA':
* read all 1700000 entries after 176s

The two lines containing fread_line and the in-line operation split have basically the same output as multiple calls to fread_string, but are slightly faster

We introduced a check for whether the string-to-number conversion works. Disabling this checks reduces the runtime from 176 seconds to 151 seconds.

We use Index0 to keep track of the number of valid entries in the tuples. We like to add 0 to index variable names as a reminder that this index is 0-based, that is the first referenced element has the index 0. In HALCON, confusion might come from the fact that iconic object indices are 1-based while tuple a vector indices are 0-based.

Try to use HALCON operators on large tuples

HALCON operators are much faster if they process multiple data in a single call instead of multiple calls that process only a single entry. Example:

dev_update_off ()
NumEntries := 100000
SomeLargeTuple1 := rand(NumEntries)
SomeLargeTuple2 := rand(NumEntries)
* slow (1-2 seconds)
for Index0 := 0 to NumEntries-1 by 1
    Sum[Index0] := SomeLargeTuple1[Index0] + SomeLargeTuple2[Index0]
endfor
* fast (3 milliseconds)
Sum := SomeLargeTuple1 + SomeLargeTuple2

In our case, the number() and regexp_replace operations are expensive if called for each point in the point cloud separately. Therefore, the script can be tuned further by collecting all coordinates and gray values in a large tuple P, and then call the number() and regexp_replace operations only once for this large tuple. Afterwards the numbers have to be deinterleaved to X, Y, Z, Gray tuples:

read_point_cloud_simple_v6.hdev

dev_update_off ()

Filename := '../tire-heindl-solutions.txt'
NumColumns := 4

count_seconds (Then)
open_file (Filename, 'input', FileHandle)

* choose LargeNum at least the size of the expected number of coordinates,
* thus, no rellocations are necessary during insertions into the tuple P.
* Expected memory allocated: LargeNum * 4 * 8 bytes = 305 MB

LargeNum := NumColumns*10000000
P := gen_tuple_const(LargeNum, 0)
Index0 := 0
while (true)
    fread_line (FileHandle, OutLine, IsEOF)
    if (IsEOF)
        break
    endif
    SplitString := split(OutLine, ' ')
    P[Index0:Index0+NumColumns-1] := SplitString
    Index0 := Index0 + NumColumns
endwhile
* remember that we initialized the tuple P probably too large, shrink-to-fit:
P := P[0:Index0-1]

* trim white-space left and right
* trim new line character
* trim leading zeros (which would indicate an octal number to 'tuple_number')
Number := number(regexp_replace(P,'^\\s*0*(.+?)\\s*\n*$', '$1'))

* check if all coordinates/gray values are reasonable:
IndexIsString := find(type_elem(Number),H_TYPE_STRING)
if (IndexIsString > -1)
    throw ('could not convert "' + Number[IndexIsString] + '" to a number')
endif

X := P[[0:NumColumns:|P|-1]]
Y := P[[1:NumColumns:|P|-1]]
Z := P[[2:NumColumns:|P|-1]]
Gray := P[[3:NumColumns:|P|-1]]

* explicit free memory (only makes sense for very large tuples)
P := HNULL

close_file (FileHandle)
count_seconds (Now)
DurationSeconds := Now - Then
Msg := 'opening file ' + Filename + ' in ' + DurationSeconds + 's'
dev_inspect_ctrl (Msg)

* runtime on machine 'MEGA':
* read all 1700000 entries after 92s

Call `split` only once

Even the split operation can be deferred to a single call after all lines are read: read_point_cloud_simple_v7.hdev

dev_update_off ()

Filename := '../tire-heindl-solutions.txt'
NumColumns := 4

count_seconds (Then)
open_file (Filename, 'input', FileHandle)

* choose LargeNum at least the size of the expected number of coordinates,
* thus, no rellocations are necessary during insertions into the tuple P.
* Expected memory allocated: LargeNum * 4 * 8 bytes = 305 MB

LargeNum := 10000000
P := gen_tuple_const(LargeNum, 0)
Index0 := 0
IsEOF := false
repeat
    fread_line (FileHandle, OutLine, IsEOF)
    P[Index0] := OutLine
    Index0 := Index0 + 1
until (IsEOF)
* remember that we initialized the tuple P probably too large, shrink-to-fit:
* we have to remove one additional entry as we added the unwanted empty last line, too
P := P[0:Index0-1-1]

P := split(P, ' ')


* trim white-space left and right
* trim new line character
* trim leading zeros (which would indicate an octal number to 'tuple_number')
Number := number(regexp_replace(P,'^\\s*0*(.+?)\\s*\n*$', '$1'))

* explicit free memory (only makes sense for very large tuples)
P := HNULL

* check if all coordinates/gray values are reasonable:
IndexIsString := find(type_elem(Number),H_TYPE_STRING)
if (IndexIsString > -1)
    throw ('could not convert "' + Number[IndexIsString] + '" to a number')
endif

X := Number[[0:NumColumns:|Number|-1]]
Y := Number[[1:NumColumns:|Number|-1]]
Z := Number[[2:NumColumns:|Number|-1]]
Gray := Number[[3:NumColumns:|Number|-1]]

* explicit free memory (only makes sense for very large tuples)
Number := HNULL


close_file (FileHandle)
count_seconds (Now)
DurationSeconds := Now - Then
Msg := 'opening file ' + Filename + ' in ' + DurationSeconds + 's'
dev_inspect_ctrl (Msg)

* runtime on machine 'MEGA':
* read all 1700000 entries after 67s

Vector operations

Using vector operations - in this case reading the lines of the file directly into a HDevelop vector - can be a further performance improvement:

read_point_cloud_simple_v8.hdev

dev_update_off ()

Filename := '../tire-heindl-solutions.txt'
NumColumns := 4

count_seconds (Then)
open_file (Filename, 'input', FileHandle)

VecOutLine.clear()
repeat
    fread_line (FileHandle, VecOutLine.at(VecOutLine.length()), IsEOF)
until (IsEOF)
convert_vector_to_tuple (VecOutLine, P)
* we have to remove one additional entry as we added the unwanted empty last line, too
P := P[0:|P|-2]
P := split(P, ' ')

* trim white-space left and right
* trim new line character
* trim leading zeros (which would indicate an octal number to 'tuple_number')
Number := number(regexp_replace(P,'^\\s*0*(.+?)\\s*\n*$', '$1'))

* explicit free memory (only makes sense for very large tuples)
P := HNULL

* check if all coordinates/gray values are reasonable:
IndexIsString := find(type_elem(Number),H_TYPE_STRING)
if (IndexIsString > -1)
    throw ('could not convert "' + Number[IndexIsString] + '" to a number')
endif

X := Number[[0:NumColumns:|Number|-1]]
Y := Number[[1:NumColumns:|Number|-1]]
Z := Number[[2:NumColumns:|Number|-1]]
Gray := Number[[3:NumColumns:|Number|-1]]

* explicit free memory (only makes sense for very large tuples)
Number := HNULL

close_file (FileHandle)
count_seconds (Now)
DurationSeconds := Now - Then
Msg := 'opening file ' + Filename + ' in ' + DurationSeconds + 's'
dev_inspect_ctrl (Msg)

* runtime on machine 'MEGA':
* read all 1700000 entries after 44s

Extract into Procedure

It is often a good idea to refactor reusable code into its own procedure. Therefore, we extract most of the code into a procedure called ReadPointCloudFromTextfile.

Then one could consider to enable JIT compiled procedure execution (in HDevelop: General Options / Experienced User / Execute procedures JIT-compiled) to get a finally very fast runtime of 21 seconds. Because of multiple issues with JIT compiling we do not generally recommend using it, at least during development.

dev_update_off ()

Filename := '../tire-heindl-solutions.txt'
NumColumns := 4

* enable JIT compiled procedure execution
* General Options / Experienced User / Execute procedures JIT-compiled

count_seconds (Then)
ReadPointCloudFromTextfile (Filename, NumColumns, X, Y, Z, Gray)
count_seconds (Now)
DurationSeconds := Now - Then
Msg := 'opening file ' + Filename + ' in ' + DurationSeconds + 's'
dev_inspect_ctrl (Msg)

* runtime on machine 'MEGA':
* read all 1700000 entries after 21s

* visualize the point cloud:
gen_object_model_3d_from_points (X, Y, Z, ObjectModel3D)
set_object_model_3d_attrib_mod (ObjectModel3D, '&gray', 'points', Gray)
X := HNULL
Y := HNULL
Z := HNULL
Gray := HNULL
dev_get_window (WindowHandle)
visualize_object_model_3d (WindowHandle, ObjectModel3D, [], [], 'intensity', '&gray', [], [], [], PoseOut)


Procedure ReadPointCloudFromTextfile ( : : Filename, NumColumns : X, Y, Z, Grayval )
VecOutLineOut.clear()

open_file (Filename, 'input', FileHandle)
repeat
    fread_line (FileHandle, VecOutLineOut.at(VecOutLineOut.length()), IsEOF)
until (IsEOF)
close_file (FileHandle)

convert_vector_to_tuple (VecOutLineOut, P)
* we would have to remove the last entry as we added the unwanted empty last line, too
* but this empty last line will be taken care of by the split() operand
P := split(P, ' ')

* trim white-space left and right
* trim new line character
* trim leading zeros (which would indicate an octal number to 'tuple_number')
Number := number(regexp_replace(P,'^\\s*0*(.+?)\\s*\n*$', '$1'))

* check if all coordinates/gray values are reasonable:
IndexIsString := find(type_elem(Number),H_TYPE_STRING)
if (IndexIsString > -1)
    throw ('could not convert "' + Number[IndexIsString] + '" to a number')
endif

X := Number[[0:NumColumns:|Number|-1]]
Y := Number[[1:NumColumns:|Number|-1]]
Z := Number[[2:NumColumns:|Number|-1]]
Gray := Number[[3:NumColumns:|Number|-1]]
return ()

Downloads

Summary

We have shown how a text-based point cloud file can be read fast, robust, and efficiently in HDevelop. Please contact us if you have any questions or need help with your project.

(Fast) parsing of a point cloud text file

2019-09-19 by Andreas Heindl