from gvsig import *
from geom import *
import commonsdialog

def main(*args):

    """Calculate Sinuosity Index and make
    reclassification following Schumm (1963) criteria"""
    #DATA:
    #input: River, roads, etc layer type LINE.
        #This input shapefile must be selected in TOC
        #The EPSG must not be projected
        
    #output: same layer with 8 new fields:
        # X0 = X coordinate field of the starting point
        # Y0 = Y coordinate field of the starting point
        # Xn = X coordinate field of the ending point
        # Yn = Y coordinate field of the ending point
        # Dist_lin = distance field between starting point and ending point
        # Dist_real = real distance field
        # Si = Sinuosity Index field
        # Schumm = reclassification field for Sinuosity Index, by Schumm

    #capa de lineas sobre la que se trabaja
    layer = currentLayer()
    features = layer.features()
    schema = layer.getSchema()

    #verify that the fields to create no longer exist in the table
    if "X0" in schema.getAttrNames() or "Y0" in schema.getAttrNames() or "Xn" in schema.getAttrNames() or "Yn" in schema.getAttrNames() or "Dist_lin" in schema.getAttrNames() or "Dist_real" in schema.getAttrNames() or "Si" in schema.getAttrNames() or "Schumm" in schema.getAttrNames():
        commonsdialog.msgbox("At least one of the fields to be created already exists:'\n'X0'\n'Y0'\n'Xn'\n'Yn'\n'Dist_lin'\n'Dist_real'\n'Si'\n'Schumm'\n'Change the name of the existing field that already exists.", "Attention",2)
    else:
            #create new schema to add necesary fields
            newSchema=createSchema(schema)
            # add X coordinate field of the starting point
            newSchema.append("X0","FLOAT",size=30, precision=10)
            # add Y coordinate field of the starting point
            newSchema.append("Y0","FLOAT",size=30, precision=10)
            # add X coordinate field of the ending point
            newSchema.append("Xn","FLOAT",size=30, precision=10)
            # add Y coordinate field of the ending point
            newSchema.append("Yn","FLOAT",size=30, precision=10)
            # add distance field between starting point and ending point
            newSchema.append("Dist_lin","FLOAT",size=30, precision=10)
            # add real distance field
            newSchema.append("Dist_real","FLOAT",size=30, precision=10)
            # add Sinuosity Index field
            newSchema.append("Si","FLOAT",size=30, precision=10)
            # add reclassification field for Sinuosity Index, by Schumm
            newSchema.append("Schumm","STRING",size=15)
            #edit layer to update schema
            layer.edit()
            layer.updateSchema(newSchema)
            layer.commit()
            
    #for each feature in the layer calculate...
    for feature in layer.features():
            #obtain first coordinate
            punto = feature.geometry().getVertex(0)
            #obtain last coordinate
            puntofin = feature.geometry().getVertex(feature.geometry().getNumVertices()- 1)
            #calculate distance between starting point and ending point
            distancia = pow((float(puntofin.getX())-float(punto.getX()))**2 + (float(puntofin.getY())-float(punto.getY()))** 2, 0.5)
            #calculate real distance
            dist_real = feature.geometry().perimeter()
            #calculate Sinuosity Index
            Ind_Sin = float(dist_real)/float(distancia)

            #edit layer and update feature values
            feature.edit()
            #obtain X coordinate of the starting point
            feature.set("X0", float(punto.getX()))
            #obtain Y coordinate of the starting point
            feature.set("Y0", float(punto.getY()))
            #obtain X coordinate of the ending point
            feature.set("Xn", float(puntofin.getX()))
            #obtain Y coordinate of the ending point
            feature.set("Yn", float(puntofin.getY()))
            #update features with calculated values: distancia
            feature.set("Dist_lin", float(distancia))
            #update features with calculated values: dist_real
            feature.set("Dist_real", float(dist_real))
            #update features with calculated values: Ind_Sin
            feature.set("Si", float (Ind_Sin))

            #reclassify the Sinuosity Index by Schumm 
            if Ind_Sin <= float(1.2) : 
                feature.set("Schumm", str("Rectilinear"))
            elif Ind_Sin > float(1.2) and Ind_Sin <= float(1.5) :
                feature.set("Schumm", str("Transitional"))
            elif Ind_Sin > float(1.5) and Ind_Sin <= float(1.7) :
                feature.set("Schumm", str("Regular"))
            elif Ind_Sin > float(1.7) and Ind_Sin <= float(2.1) :
                feature.set("Schumm", str("Irregular"))
            elif Ind_Sin > float(2.1) :
                feature.set("Schumm", str("Tortuous"))

            #update and save data edited
            layer.update(feature)
            layer.commit()