<jupyter_start><jupyter_text>Select points<jupyter_code>import numpy as np
import matplotlib.pyplot as plt
import osgeo.gdal as gdal
import os
from maxvol_cut import rect_maxvol_cut, f_no_cut, f_penal_2D
from tools import norm_data, add_coords, gen_input, extend_score, points_selection, f_no_cut, f_cut_eps, calc_score, good_points_brute_force, idx_to_idx

%matplotlib inline

# Set BIG figure
plt.rcParams["figure.figsize"] = (12, 8)
dr = '../local_tifs'  
NUMBER_OF_POINTS=16
fl_names = list(filter(lambda fl: fl.endswith('.tif'), os.listdir(dr)))
files = list(map(lambda x: gdal.Open(os.path.join(dr, x)), fl_names[::-1]))
arrays = list(map(lambda x: x.ReadAsArray().flatten(), files))
shapes = [x.ReadAsArray().shape for x in files]
nodatas = list(map(lambda x: x.GetRasterBand(1).GetNoDataValue(), files))
names = list(map(lambda x: x.replace('.tif','').split('.')[0], fl_names))
dem_raw = gdal.Open('../10_3857.tif')
dem = dem_raw.ReadAsArray()
dem_flat = dem.flatten()
dem_nodata = dem_raw.GetRasterBand(1).GetNoDataValue()
init_dem_shape = dem.shape
#delete nodata
idx_nodata_0 = np.where(arrays[0] == nodatas[0])[0]
arrays_no_nodatas = np.zeros((len(arrays[0])-len(idx_nodata_0), len(arrays)))

idx_dem_nodata = np.where(dem_flat == dem_nodata)[0]
idx_dem = np.where(arrays[0] != nodatas[0])[0]
dem_no_nodata = np.delete(dem_flat, idx_nodata_0)

for i in range(len(arrays)):
    idx_nodata = np.where(arrays[i] == nodatas[i])[0]
    array = arrays[i].copy()
    array[idx_nodata]=0
    arrays_no_nodatas[:,i]  = np.delete(array, idx_nodata_0)
data_arr = arrays_no_nodatas.copy()
# Prepare data
# U can normilize data, and/or add coords to it

mode = 1 # Change to 0, 1, 2 or 3
X, fn_X_embedded = gen_input(mode, data_arr, shapes, idx_dem)
X.shape
#function for distance between points
f_cut = lambda idx, i : f_cut_eps(idx, i, X=X, eps=0.3)

#function for distence from border
f_penal = f_penal_2D(X = X[:, -2], Y = X[:, -1], bnd = 0.3, level = 0.3)
#result of selection
#if result print "Failed", the number of chosing points will be smaller than NUMBER_OF_POINTS
#to fix it, reduce eps in f_cut and bnd, level in f_penal

result = points_selection(X, n_pnts = NUMBER_OF_POINTS, cut_fun = f_cut, penalty = None)<jupyter_output><empty_output><jupyter_text>Export coordinates To export coordinates, you should know border values of your DEM (xmin, ymin, xmax, ymax).They could be copied from QGIS, features of layer<jupyter_code>xmin, ymin, xmax, ymax = [37.7928399298814384,51.9023655604198453,37.8064010466814366,51.9077426838598441]

dem_flat_img = dem_flat.copy()-np.min(dem_flat)
dem_flat_img[np.where(dem_flat == dem_nodata)] = float('NaN')
st = dem_flat_img.reshape(init_dem_shape)

xv = np.linspace(xmin,xmax, num=st.shape[1])
yv = np.linspace(ymax,ymin, num=st.shape[0])
coords = np.meshgrid(xv,yv)

mask = idx_dem


#select corresponding points by indecies
y_c,x_c = coords[0].flatten()[mask, None],coords[1].flatten()[mask, None]
y_idx, x_idx = y_c[result],x_c[result]
coord_idx = np.hstack((y_idx,x_idx))

np.savetxt('maxvol_result_'+str(NUMBER_OF_POINTS)+'.csv', coord_idx, delimiter=';')<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1445859872611465, 'fraction_numerical': 0.03821656050955414, 'mean_word_length': 6.270833333333333, 'pattern_counts': {'":': 0, '<': 9, '<?xml version=': 0, '>': 9, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '2', 'dir': 'jupyter-structured-clean-dedup', 'id': '42403', 'license': 'MIT', 'n_tokens_mistral': 1338, 'n_tokens_neox': 1264, 'n_words': 274, 'path': 'jupyter_notebooks/maxvol_points_selection.ipynb', 'repo_name': 'EDSEL-skoltech/maxvol_sampling', 'size': '3140'}