Skip to content

Commit 68a271f

Browse files
debpaldebpal
committed
code upgrade
1 parent 3a96d8a commit 68a271f

File tree

4 files changed

+750
-232
lines changed

4 files changed

+750
-232
lines changed

Diff for: ‎GeoAnalyze/raster.py

+286-28
Original file line numberDiff line numberDiff line change
@@ -1,6 +1,10 @@
11
import rasterio
2+
import rasterio.features
3+
import rasterio.mask
4+
import geopandas
25
import pandas
36
import numpy
7+
import typing
48

59

610
class Raster:
@@ -9,7 +13,7 @@ class Raster:
913
Provides functionality for raster file operations.
1014
'''
1115

12-
def count_non_nodata_cells(
16+
def count_data_cells(
1317
self,
1418
input_file: str
1519
) -> int:
@@ -61,50 +65,110 @@ def count_nodata_cells(
6165

6266
return output
6367

64-
def percentage_unique_integers(
68+
def counting_unique_values(
6569
self,
66-
input_file: str
70+
input_file: str,
71+
csv_file: str,
72+
multiplier: float = 1
6773
) -> pandas.DataFrame:
6874

6975
'''
70-
Computes the percentage of unique integer values in the raster array.
76+
Computes the count of unique data values in a raster array. If the data contains decimal values,
77+
the input multiplier is used to scale the decimal values to integers for counting purposes.
78+
The data are then scaled back to the decimal values by dividing by the multiplier.
7179
7280
Parameters
7381
----------
7482
input_file : str
7583
Path to the input raster file.
7684
85+
csv_file : str
86+
Path to save the output csv file.
87+
88+
multiplier : float, optional
89+
A factor to multiply raster values to handle decimal values by rounding.
90+
Default is 1, which implies no scaling.
91+
7792
Returns
7893
-------
7994
DataFrame
80-
A DataFrame containing the unique integer values in the raster,
81-
their counts, and their counts as a percentage of the total.
95+
A DataFrame containing the raster values, their counts,
96+
and their counts as a percentage of the total.
8297
'''
8398

8499
with rasterio.open(input_file) as input_raster:
85100
raster_profile = input_raster.profile
86-
if 'int' not in raster_profile['dtype']:
87-
raise Exception('Input raster data must be integer type.')
88-
else:
89-
# DataFrame
90-
raster_array = input_raster.read(1)
91-
valid_array = raster_array[raster_array != raster_profile['nodata']]
92-
value, count = numpy.unique(
93-
valid_array,
94-
return_counts=True
95-
)
96-
df = pandas.DataFrame({'Value': value, 'Count': count})
97-
df['Percent(%)'] = 100 * df['Count'] / df['Count'].sum()
101+
raster_array = input_raster.read(1)
102+
value_array = (multiplier * raster_array[raster_array != raster_profile['nodata']]).round()
103+
value, count = numpy.unique(
104+
value_array,
105+
return_counts=True
106+
)
107+
df = pandas.DataFrame({'Value': value, 'Count': count})
108+
df['Value'] = df['Value'] / multiplier
109+
df['Count(%)'] = 100 * df['Count'] / df['Count'].sum()
110+
df['Cumulative_Count(%)'] = df['Count(%)'].cumsum()
111+
df.to_csv(
112+
path_or_buf=csv_file,
113+
index_label='Index',
114+
float_format='%.2f'
115+
)
98116

99117
return df
100118

119+
def boundary_polygon(
120+
self,
121+
input_file: str,
122+
output_file: str
123+
) -> geopandas.GeoDataFrame:
124+
125+
'''
126+
Extracts boundary polygons from a raster array.
127+
128+
Parameters
129+
----------
130+
input_file : str
131+
Path to the input raster file.
132+
133+
output_file : str
134+
Path to save the output shapefile.
135+
136+
Returns
137+
-------
138+
GeoDataFrame
139+
A GeoDataFrame containing the boundary polygons extracted from the raster.
140+
'''
141+
142+
with rasterio.open(input_file) as input_raster:
143+
raster_array = input_raster.read(1)
144+
raster_array[raster_array != input_raster.nodata] = 1
145+
mask = raster_array == 1
146+
boundary_shapes = rasterio.features.shapes(
147+
source=raster_array,
148+
mask=mask,
149+
transform=input_raster.transform
150+
)
151+
boundary_features = [
152+
{'geometry': geom, 'properties': {'value': val}} for geom, val in boundary_shapes
153+
]
154+
gdf = geopandas.GeoDataFrame.from_features(
155+
features=boundary_features,
156+
crs=input_raster.crs
157+
)
158+
gdf = gdf[gdf.is_valid].reset_index(drop=True)
159+
gdf['bid'] = range(1, gdf.shape[0] + 1)
160+
gdf = gdf[['bid', 'geometry']]
161+
gdf.to_file(output_file)
162+
163+
return gdf
164+
101165
def rescaling_resolution(
102166
self,
103167
input_file: str,
104168
target_resolution: int,
105169
resampling_method: str,
106170
output_file: str
107-
) -> float:
171+
) -> list[list[float]]:
108172

109173
'''
110174
Rescales the raster array from the existing resolution to a new target resolution.
@@ -128,8 +192,9 @@ def rescaling_resolution(
128192
129193
Returns
130194
-------
131-
float
132-
The resolution of the rescaled raster array.
195+
list
196+
A nested list representation of the 3x3 affine transformation matrix
197+
of the reprojected raster array. Each sublist represents a row of the matrix.
133198
'''
134199

135200
# mapping of resampling methods
@@ -179,17 +244,19 @@ def rescaling_resolution(
179244
dst_crs=input_raster.crs,
180245
resampling=resampling_function[resampling_method]
181246
)
182-
output_resolution = float(rescaled_raster[1][0])
247+
affine_matrix = [
248+
list(rescaled_raster[1])[i:i + 3] for i in [0, 3, 6]
249+
]
183250

184-
return output_resolution
251+
return affine_matrix
185252

186253
def rescaling_resolution_with_mask(
187254
self,
188255
input_file: str,
189256
mask_file: str,
190257
resampling_method: str,
191258
output_file: str
192-
) -> float:
259+
) -> list[list[float]]:
193260

194261
'''
195262
Rescales the raster array from its existing resolution
@@ -214,8 +281,9 @@ def rescaling_resolution_with_mask(
214281
215282
Returns
216283
-------
217-
float
218-
The resolution of the rescaled raster array.
284+
list
285+
A nested list representation of the 3x3 affine transformation matrix
286+
of the reprojected raster array. Each sublist represents a row of the matrix.
219287
'''
220288

221289
# mapping of resampling methods
@@ -270,6 +338,196 @@ def rescaling_resolution_with_mask(
270338
dst_crs=mask_raster.crs,
271339
resampling=resampling_function[resampling_method]
272340
)
273-
output_resolution = float(rescaled_raster[1][0])
341+
affine_matrix = [
342+
list(rescaled_raster[1])[i:i + 3] for i in [0, 3, 6]
343+
]
344+
345+
return affine_matrix
346+
347+
def reproject_crs(
348+
self,
349+
input_file: str,
350+
resampling_method: str,
351+
target_crs: str,
352+
output_file: str
353+
) -> list[list[float]]:
354+
355+
'''
356+
Reprojects a raster array to a new Coordinate Reference System.
357+
358+
Parameters
359+
----------
360+
input_file : str
361+
Path to the input raster file.
362+
363+
resampling_method : str
364+
Method used for raster resampling. Supported options are:
365+
- 'nearest': Nearest-neighbor interpolation.
366+
- 'bilinear': Bilinear interpolation.
367+
- 'cubic': Cubic convolution interpolation.
368+
369+
target_crs : str
370+
Target Coordinate Reference System for the output raster (e.g., 'EPSG:4326').
371+
372+
output_file : str
373+
Path to save the reprojected raster file.
374+
375+
Returns
376+
-------
377+
list
378+
A nested list representation of the 3x3 affine transformation matrix
379+
of the reprojected raster array. Each sublist represents a row of the matrix.
380+
'''
381+
382+
# mapping of resampling methods
383+
resampling_function = {
384+
'nearest': rasterio.enums.Resampling.nearest,
385+
'bilinear': rasterio.enums.Resampling.bilinear,
386+
'cubic': rasterio.enums.Resampling.cubic
387+
}
388+
389+
# check resampling method
390+
if resampling_method in resampling_function.keys():
391+
pass
392+
else:
393+
raise Exception('Input resampling method is not supported.')
394+
395+
# reproject Coordinate Reference System
396+
with rasterio.open(input_file) as input_raster:
397+
raster_profile = input_raster.profile
398+
# output raster parameters
399+
output_transform, output_width, output_height = rasterio.warp.calculate_default_transform(
400+
src_crs=input_raster.crs,
401+
dst_crs=target_crs,
402+
width=input_raster.width,
403+
height=input_raster.height,
404+
left=input_raster.bounds.left,
405+
bottom=input_raster.bounds.bottom,
406+
right=input_raster.bounds.right,
407+
top=input_raster.bounds.top
408+
)
409+
# output raster profile
410+
raster_profile.update(
411+
{
412+
'transform': output_transform,
413+
'width': output_width,
414+
'height': output_height,
415+
'crs': target_crs
416+
}
417+
)
418+
# saving output raster
419+
with rasterio.open(output_file, 'w', **raster_profile) as output_raster:
420+
reproject_raster = rasterio.warp.reproject(
421+
source=rasterio.band(input_raster, 1),
422+
destination=rasterio.band(output_raster, 1),
423+
src_transform=input_raster.transform,
424+
src_crs=input_raster.crs,
425+
dst_transform=output_transform,
426+
dst_crs=target_crs,
427+
resampling=resampling_function[resampling_method]
428+
)
429+
affine_matrix = [
430+
list(reproject_raster[1])[i:i + 3] for i in [0, 3, 6]
431+
]
432+
433+
return affine_matrix
434+
435+
def set_value_to_nodata(
436+
self,
437+
input_file: str,
438+
target_value: list[float],
439+
nodata: float,
440+
output_file: str
441+
) -> dict[str, typing.Any]:
442+
443+
'''
444+
Converts specified values in a raster array to nodata.
445+
446+
Parameters
447+
----------
448+
input_file : str
449+
Path to the input raster file.
450+
451+
target_value : list
452+
List of values in the input raster array to convert to nodata.
453+
454+
nodata : float
455+
The nodata value to assign in the output raster.
456+
457+
output_file : str
458+
Path to save the output raster file.
459+
460+
Returns
461+
-------
462+
dict
463+
A profile dictionary containing metadata about the output raster.
464+
'''
465+
466+
with rasterio.open(input_file) as input_raster:
467+
raster_profile = input_raster.profile
468+
input_array = input_raster.read(1)
469+
output_array = numpy.where(
470+
numpy.isin(input_array, target_value),
471+
nodata,
472+
input_array
473+
)
474+
raster_profile.update(
475+
{
476+
'nodata': nodata
477+
}
478+
)
479+
with rasterio.open(output_file, 'w', **raster_profile) as output_raster:
480+
output_raster.write(output_array, 1)
481+
output_profile = dict(output_raster.profile)
482+
483+
return output_profile
484+
485+
def clipping_by_shapes(
486+
self,
487+
input_file: str,
488+
shape_file: str,
489+
output_file: str
490+
) -> dict[str, typing.Any]:
491+
492+
'''
493+
Clips a raster file using a given shape file.
494+
495+
Parameters
496+
----------
497+
input_file : str
498+
Path to the input raster file.
499+
500+
shape_file : str
501+
Path to the input shape file used for clipping.
502+
503+
output_file : str
504+
Path to save the output raster file.
505+
506+
Returns
507+
-------
508+
dict
509+
A profile dictionary containing metadata about the output raster.
510+
'''
511+
512+
with rasterio.open(input_file) as input_raster:
513+
raster_profile = input_raster.profile.copy()
514+
gdf = geopandas.read_file(shape_file)
515+
gdf = gdf.to_crs(str(raster_profile['crs']))
516+
output_array, output_transform = rasterio.mask.mask(
517+
dataset=input_raster,
518+
shapes=gdf.geometry.tolist(),
519+
all_touched=True,
520+
crop=True
521+
)
522+
raster_profile.update(
523+
{
524+
'height': output_array.shape[1],
525+
'width': output_array.shape[2],
526+
'transform': output_transform
527+
}
528+
)
529+
with rasterio.open(output_file, 'w', **raster_profile) as output_raster:
530+
output_raster.write(output_array)
531+
output_profile = dict(output_raster.profile)
274532

275-
return output_resolution
533+
return output_profile

0 commit comments

Comments
 (0)