main.py

from time import time
import matplotlib.path as mplPath
import numpy as np
from progress import printDate, printProgressBar, printTiming

# this file contains boundaries for AIN TEMOUCHENTE
BOUNDS_FILE = 'bounds.txt'
# this two files are the patches covering the AIN TEMOUCHENTE area
XYZ_PATCHS = ['coordinates/N35W002.xyz', 'coordinates/N35W001.xyz']
# this is an output file, csv format so it could be used in mapinfo
OUTPUT_DIR = 'outputs/'
OUTPUT_FILE = 'output.csv'


def readBounds(file_name):
    """
    this function loads all the coordinates and creates a polygone
    it represents the boundaries
    """
    bounds = []
    with open(file_name, 'r') as f:
        for line in f:
            coordinates = line.strip().split(' ')
            bounds.append([float(coordinates[0]), float(coordinates[1])])

    return mplPath.Path(np.array(bounds))


def extractXYZ(polygon):
    """    
    this function reads a provided patch and extracts the ones that
    are contained inside the polygon created by the readBounds function
    """
    TOTAL_LINES = 2880000
    currentLine = 0
    for idx, patch in enumerate(XYZ_PATCHS):
        mode = 'w' if idx == 0 else'a'
        input = open(patch, "r")
        output = open(OUTPUT_DIR+OUTPUT_FILE, mode)
        # print file headers
        if(mode is'w'):
            print("\"longitude\",\"latitude\",\"altitude\"", file=output, end='\n')
        for line in input:
            coordinates = line.strip().split(' ')
            coordinates = (float(coordinates[0]), float(coordinates[1]))
            printProgressBar(currentLine, TOTAL_LINES)
            currentLine += 1
            if(polygon.contains_point(coordinates)):
                line = ','.join(map(str, line.strip().split(' ')))
                print(line, file=output, end='\n')
        input.close()
        output.close()


def extractLinesEvery(lineNumber):
    """
    this function is to extract a line every specifique number
    from the file generated by the extractXYZ and outputs the same CSV
    """
    index = 0
    with open(OUTPUT_DIR+OUTPUT_FILE, 'r') as input:
        points = []
        with open(OUTPUT_DIR+str(lineNumber)+'%_'+OUTPUT_FILE, 'w') as output:
            for line in input:
                if(index % lineNumber == 0):
                    print(line, file=output, end='')


if __name__ == '__main__':
    start = time()
    printDate('Extraction START TIME: ')
    polygon = readBounds(BOUNDS_FILE)
    extractXYZ(polygon)
    printDate('\nExtraction END TIME: ')
    printTiming(start)
    extractLinesEvery(10)  # lines