|
1 | | -""" |
2 | | -cBLUE (comprehensive Bathymetric Lidar Uncertainty Estimator) |
3 | | -Copyright (C) 2019 |
4 | | -Oregon State University (OSU) |
5 | | -Center for Coastal and Ocean Mapping/Joint Hydrographic Center, University of New Hampshire (CCOM/JHC, UNH) |
6 | | -NOAA Remote Sensing Division (NOAA RSD) |
7 | | -
|
8 | | -This library is free software; you can redistribute it and/or |
9 | | -modify it under the terms of the GNU Lesser General Public |
10 | | -License as published by the Free Software Foundation; either |
11 | | -version 2.1 of the License, or (at your option) any later version. |
12 | | -
|
13 | | -This library is distributed in the hope that it will be useful, |
14 | | -but WITHOUT ANY WARRANTY; without even the implied warranty of |
15 | | -MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
16 | | -Lesser General Public License for more details. |
17 | | -
|
18 | | -You should have received a copy of the GNU Lesser General Public |
19 | | -License along with this library; if not, write to the Free Software |
20 | | -Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA |
21 | | -
|
22 | | -Contact: |
23 | | -Christopher Parrish, PhD |
24 | | -School of Construction and Civil Engineering |
25 | | -101 Kearney Hall |
26 | | -Oregon State University |
27 | | -Corvallis, OR 97331 |
28 | | -(541) 737-5688 |
29 | | - |
30 | | -
|
31 | | -Last Edited By: |
32 | | -Keana Kief (OSU) |
33 | | -July 25th, 2023 |
34 | | -
|
35 | | -
|
36 | | -""" |
37 | | -import os |
38 | | -import pandas as pd |
39 | | -import logging |
40 | | -import numpy as np |
41 | | -import numexpr as ne |
42 | | -import laspy |
43 | | - |
44 | | -logger = logging.getLogger(__name__) |
45 | | - |
46 | | -""" |
47 | | -This class provides the functionality to load las files into cBLUE. One Las object |
48 | | -is created for each loaded las file. |
49 | | -""" |
50 | | - |
51 | | - |
52 | | -class Las: |
53 | | - def __init__(self, las): |
54 | | - self.las = las |
55 | | - self.las_short_name = os.path.split(las)[-1] |
56 | | - self.las_base_name = self.las_short_name.replace(".las", "") |
57 | | - self.inFile = laspy.read(self.las) |
58 | | - self.points_to_process = self.inFile.points |
59 | | - self.unq_flight_lines = self.get_flight_line_ids() |
60 | | - self.num_file_points = self.points_to_process.array.shape[0] |
61 | | - |
62 | | - """index list that would sort gps_time (to be used to |
63 | | - later when exporting las data and calculated tpu to a las |
64 | | - file |
65 | | - """ |
66 | | - self.t_argsort = None |
67 | | - |
68 | | - def get_bathy_points(self): |
69 | | - class_codes = {"BATHYMETRY": 26} |
70 | | - bathy_inds = self.inFile.raw_classification == class_codes["BATHYMETRY"] |
71 | | - return self.inFile.points.array[bathy_inds]["point"] |
72 | | - |
73 | | - def get_flight_line_ids(self): |
74 | | - """generates a list of unique flight line ids |
75 | | -
|
76 | | - This method returns a list of unique flight line ids. |
77 | | -
|
78 | | - :return: ndarray |
79 | | - """ |
80 | | - |
81 | | - # pandas' unique faster than numpy's ? |
82 | | - return pd.unique(self.points_to_process["pt_src_id"]) |
83 | | - |
84 | | - def get_flight_line(self, sensor_type): |
85 | | - """retrieves the x, y, z, and timestamp data from the las data points |
86 | | -
|
87 | | - The x, y, and z values in the las file are stored as integers. The |
88 | | - scale and offset values in the las file header are used to convert |
89 | | - the integer values to decimal values with centimeter precision. |
90 | | -
|
91 | | - :param ? fl: |
92 | | - :return: np.array, np.array, np.array, np.array |
93 | | - """ |
94 | | - |
95 | | - #xyz_to_coordinate converts the x, y, z integer values to coordinate values |
96 | | - x, y, z = self.xyz_to_coordinate() |
97 | | - |
98 | | - t = self.points_to_process["gps_time"] |
99 | | - |
100 | | - if "classification" in self.points_to_process.array.dtype.names: |
101 | | - c = self.points_to_process["classification"] |
102 | | - elif "classification_flags" in self.points_to_process.array.dtype.names: |
103 | | - c = self.points_to_process["classification_flags"] |
104 | | - else: |
105 | | - raise Exception("Unknown las version or missing classification attribute.") |
106 | | - |
107 | | - self.t_argsort = t.argsort() |
108 | | - |
109 | | - # Check if this is a multi beam sensor |
110 | | - if(sensor_type == "multi"): |
111 | | - #Get the fan angle and multiply it by 0.006 to convert to degrees |
112 | | - fan_angle = self.inFile.scan_angle*0.006 |
113 | | - #Add xyztcf to an array together |
114 | | - xyztcf = np.vstack([x, y, z, t, c, fan_angle]).T |
115 | | - # logger.las(f"{sensor_type} Fan Angle: {fan_angle}") |
116 | | - else: |
117 | | - #Fan angle is not used by the other sensors |
118 | | - #Add xyztc to an array together |
119 | | - xyztcf = np.vstack([x, y, z, t, c]).T |
120 | | - |
121 | | - flight_lines = self.points_to_process["pt_src_id"] |
122 | | - |
123 | | - return xyztcf, self.t_argsort, flight_lines |
124 | | - |
125 | | - |
126 | | - def xyz_to_coordinate(self): |
127 | | - """The x, y, and z values in the las file are stored as integers. The |
128 | | - scale and offset values in the las file header are used to convert |
129 | | - the integer values to decimal values with centimeter precision.""" |
130 | | - |
131 | | - scale_x = np.asarray(self.inFile.header.scales[0]) |
132 | | - scale_y = np.asarray(self.inFile.header.scales[1]) |
133 | | - scale_z = np.asarray(self.inFile.header.scales[2]) |
134 | | - |
135 | | - offset_x = np.asarray(self.inFile.header.offsets[0]) |
136 | | - offset_y = np.asarray(self.inFile.header.offsets[1]) |
137 | | - offset_z = np.asarray(self.inFile.header.offsets[2]) |
138 | | - |
139 | | - X = self.points_to_process["X"] |
140 | | - Y = self.points_to_process["Y"] |
141 | | - Z = self.points_to_process["Z"] |
142 | | - |
143 | | - |
144 | | - x = ne.evaluate("X * scale_x + offset_x") |
145 | | - y = ne.evaluate("Y * scale_y + offset_y") |
146 | | - z = ne.evaluate("Z * scale_z + offset_z") |
147 | | - |
148 | | - return x, y, z |
| 1 | +""" |
| 2 | +cBLUE (comprehensive Bathymetric Lidar Uncertainty Estimator) |
| 3 | +Copyright (C) 2019 |
| 4 | +Oregon State University (OSU) |
| 5 | +Center for Coastal and Ocean Mapping/Joint Hydrographic Center, University of New Hampshire (CCOM/JHC, UNH) |
| 6 | +NOAA Remote Sensing Division (NOAA RSD) |
| 7 | +
|
| 8 | +This library is free software; you can redistribute it and/or |
| 9 | +modify it under the terms of the GNU Lesser General Public |
| 10 | +License as published by the Free Software Foundation; either |
| 11 | +version 2.1 of the License, or (at your option) any later version. |
| 12 | +
|
| 13 | +This library is distributed in the hope that it will be useful, |
| 14 | +but WITHOUT ANY WARRANTY; without even the implied warranty of |
| 15 | +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
| 16 | +Lesser General Public License for more details. |
| 17 | +
|
| 18 | +You should have received a copy of the GNU Lesser General Public |
| 19 | +License along with this library; if not, write to the Free Software |
| 20 | +Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA |
| 21 | +
|
| 22 | +Contact: |
| 23 | +Christopher Parrish, PhD |
| 24 | +School of Construction and Civil Engineering |
| 25 | +101 Kearney Hall |
| 26 | +Oregon State University |
| 27 | +Corvallis, OR 97331 |
| 28 | +(541) 737-5688 |
| 29 | + |
| 30 | +
|
| 31 | +Last Edited By: |
| 32 | +Keana Kief (OSU) |
| 33 | +September 20th, 2024 |
| 34 | +
|
| 35 | +
|
| 36 | +""" |
| 37 | +import os |
| 38 | +import pandas as pd |
| 39 | +import logging |
| 40 | +import numpy as np |
| 41 | +import numexpr as ne |
| 42 | +import laspy |
| 43 | + |
| 44 | +logger = logging.getLogger(__name__) |
| 45 | + |
| 46 | +""" |
| 47 | +This class provides the functionality to load las files into cBLUE. One Las object |
| 48 | +is created for each loaded las file. |
| 49 | +""" |
| 50 | + |
| 51 | + |
| 52 | +class Las: |
| 53 | + def __init__(self, las): |
| 54 | + self.las = las |
| 55 | + self.las_short_name = os.path.split(las)[-1] |
| 56 | + if ".las" in self.las_short_name: |
| 57 | + self.las_base_name = self.las_short_name.replace(".las", "") |
| 58 | + else: |
| 59 | + self.las_base_name = self.las_short_name.replace(".laz", "") |
| 60 | + self.inFile = laspy.read(self.las) |
| 61 | + self.points_to_process = self.inFile.points |
| 62 | + self.unq_flight_lines = self.get_flight_line_ids() |
| 63 | + self.num_file_points = self.points_to_process.array.shape[0] |
| 64 | + |
| 65 | + """index list that would sort gps_time (to be used to |
| 66 | + later when exporting las data and calculated tpu to a las |
| 67 | + file |
| 68 | + """ |
| 69 | + self.t_argsort = None |
| 70 | + |
| 71 | + def get_bathy_points(self): |
| 72 | + class_codes = {"BATHYMETRY": 26} |
| 73 | + bathy_inds = self.inFile.raw_classification == class_codes["BATHYMETRY"] |
| 74 | + return self.inFile.points.array[bathy_inds]["point"] |
| 75 | + |
| 76 | + def get_flight_line_ids(self): |
| 77 | + """generates a list of unique flight line ids |
| 78 | +
|
| 79 | + This method returns a list of unique flight line ids. |
| 80 | +
|
| 81 | + :return: ndarray |
| 82 | + """ |
| 83 | + |
| 84 | + # pandas' unique faster than numpy's ? |
| 85 | + return pd.unique(self.points_to_process["pt_src_id"]) |
| 86 | + |
| 87 | + def get_flight_line(self, sensor_type): |
| 88 | + """retrieves the x, y, z, and timestamp data from the las data points |
| 89 | +
|
| 90 | + The x, y, and z values in the las file are stored as integers. The |
| 91 | + scale and offset values in the las file header are used to convert |
| 92 | + the integer values to decimal values with centimeter precision. |
| 93 | +
|
| 94 | + :param ? fl: |
| 95 | + :return: np.array, np.array, np.array, np.array |
| 96 | + """ |
| 97 | + |
| 98 | + #xyz_to_coordinate converts the x, y, z integer values to coordinate values |
| 99 | + x, y, z = self.xyz_to_coordinate() |
| 100 | + |
| 101 | + t = self.points_to_process["gps_time"] |
| 102 | + |
| 103 | + if "classification" in self.points_to_process.array.dtype.names: |
| 104 | + c = self.points_to_process["classification"] |
| 105 | + elif "classification_flags" in self.points_to_process.array.dtype.names: |
| 106 | + c = self.points_to_process["classification_flags"] |
| 107 | + else: |
| 108 | + raise Exception("Unknown las version or missing classification attribute.") |
| 109 | + |
| 110 | + self.t_argsort = t.argsort() |
| 111 | + |
| 112 | + # Check if this is a multi beam sensor |
| 113 | + if(sensor_type == "multi"): |
| 114 | + #Get the fan angle and multiply it by 0.006 to convert to degrees |
| 115 | + fan_angle = self.inFile.scan_angle*0.006 |
| 116 | + #Add xyztcf to an array together |
| 117 | + xyztcf = np.vstack([x, y, z, t, c, fan_angle]).T |
| 118 | + # logger.las(f"{sensor_type} Fan Angle: {fan_angle}") |
| 119 | + else: |
| 120 | + #Fan angle is not used by the other sensors |
| 121 | + #Add xyztc to an array together |
| 122 | + xyztcf = np.vstack([x, y, z, t, c]).T |
| 123 | + |
| 124 | + flight_lines = self.points_to_process["pt_src_id"] |
| 125 | + |
| 126 | + return xyztcf, self.t_argsort, flight_lines |
| 127 | + |
| 128 | + |
| 129 | + def xyz_to_coordinate(self): |
| 130 | + """The x, y, and z values in the las file are stored as integers. The |
| 131 | + scale and offset values in the las file header are used to convert |
| 132 | + the integer values to decimal values with centimeter precision.""" |
| 133 | + |
| 134 | + scale_x = np.asarray(self.inFile.header.scales[0]) |
| 135 | + scale_y = np.asarray(self.inFile.header.scales[1]) |
| 136 | + scale_z = np.asarray(self.inFile.header.scales[2]) |
| 137 | + # print(f"\nLAS File: {self.las_short_name}") |
| 138 | + # print(f"\nScale X: {scale_x}") |
| 139 | + # print(f"Scale Y: {scale_y}") |
| 140 | + # print(f"Scale Z: {scale_z}") |
| 141 | + |
| 142 | + offset_x = np.asarray(self.inFile.header.offsets[0]) |
| 143 | + offset_y = np.asarray(self.inFile.header.offsets[1]) |
| 144 | + offset_z = np.asarray(self.inFile.header.offsets[2]) |
| 145 | + |
| 146 | + # print(f"\nOffset X: {offset_x}") |
| 147 | + # print(f"Offset Y: {offset_y}") |
| 148 | + # print(f"Offset Z: {offset_z}") |
| 149 | + |
| 150 | + |
| 151 | + X = self.points_to_process["X"] |
| 152 | + Y = self.points_to_process["Y"] |
| 153 | + Z = self.points_to_process["Z"] |
| 154 | + |
| 155 | + # print(f"\nX: {X}") |
| 156 | + # print(f"Y: {Y}") |
| 157 | + # print(f"Z: {Z}") |
| 158 | + |
| 159 | + |
| 160 | + x = ne.evaluate("X * scale_x + offset_x") |
| 161 | + y = ne.evaluate("Y * scale_y + offset_y") |
| 162 | + z = ne.evaluate("Z * scale_z + offset_z") |
| 163 | + |
| 164 | + # print(f"\nx: {x}") |
| 165 | + # print(f"y: {y}") |
| 166 | + # print(f"z: {z}") |
| 167 | + |
| 168 | + return x, y, z |
0 commit comments