

from os.path import basename
#from more_itertools import unique_everseen #make sure more_itertools is installed (conda install more-itertools)
from math import sqrt
from math import floor



import argparse


parser = argparse.ArgumentParser()
requiredNamed = parser.add_argument_group('required arguments')
#requiredNamed.add_argument("-d", dest='threshold', help="distance threshold for atoms near waters", type=float, default=4.0, required=True)
requiredNamed.add_argument("-i", dest='structure', help="input file, expected to be in .pdb format", type=str, required=True)
parser.add_argument("-e", dest='Even', help="boolean to set the dimensions to even by default", default=False, type=bool)
args = parser.parse_args()

filenames = []


filenames.append(basename(args.structure[:-4]))

print(filenames)
#print args.threshold
#structure = sys.argv[1]
#distance_cutoff = sys.argv[2]
# HOLDS IN INFORMATION OF ATOM
ATOM = []
# HOLDS IN ALL NON HOH ATOMS THAT ARE LABELED HETAM
HETAM = []
# HOLDS ALL THE HOH ATOMS
HETAMHOH = []
DIS = []
Unique = []
newHOH =[]  
Vectors = []
Dictionary = {}




#POINT CLASS USED TO HOLD POINTS OF ALL DIFFERENT ATOMS
class Point:
    #STORES THE X Y AND Z COORDINATES
        def __init__(Self,x, y , z ):
                Self.x = x
                Self.y = y
                Self.z = z
        #FINDS DISTANCE BETWEEN TWO POINTS
        def Distance(firstpoint, secondPoint):
        	d = sqrt(((firstpoint.x-secondPoint.x)**2)+
                    ((firstpoint.y-secondPoint.y)**2)+
                    ((firstpoint.z-secondPoint.z)**2))
        	return d
        def length(Self):
 			return sqrt((float(Self.x)**2)+(float(Self.y)**2)+(float(Self.z)**2))

		#SHOWS OUTPUT
        def __str__(Self):
                return(str(Self.x) + "," + str(Self.y)+ "," + str(Self.z))

class Atom:
  #ATOM CLASS HOLDS INFORMATION READ IN FROM FILE

  def __init__(Self,name,idnum,char,tp,att1,att2,x,y,z,av,bv,n,line):
    Self.name = name
    Self.idnum = idnum
    Self.Element = char
    Self.res = tp
    Self.chain = att1
    Self.resnum = att2
    Self.point = Point(float(x),float(y),float(z))
    Self.AValue = av
    Self.BValue = bv
    Self.AtomN = n
    Self.line = line
    Self.Acceptor = True
    Self.HOHDIS = ""
    Self.Associate = ""
    Self.nearHeavyAtom = ""
    Self.nearHydrogen = []
  def setnearHydrogen(Self, H):
    Self.nearHydrogen.append(H)
  def setnearHeavyAtom(Self, HA):
    Self.nearHeavyAtom = HA
  def setHOHDIS(Self, HOHnum):
    Self.HOHDIS = HOHnum
  def setAssociate(Self, HOH):
    Self.Associate = HOH
  def getPoint(Self):
    return Self.point
  def changeAcceptor(Self):
    Self.Acceptor = False
  def getLine(Self):
    return Self.line
  def checkSame(self, Atom):
    # Simple check function that will return true if two atoms have the same name and same coordinates
    if self.name == Atom.name and self.point == Atom.point:
        return True
    else: return False
  def __str__(Self):
    return Self.line

def parseData(line):
    if "CRYST1" in line or "REMARK" in line or "END" in line: 
        print("Skipping comment line")
        return
    if not('CL') in line and not('NA') in line and not('MG') in line and not('SO4') in line and not('PO4') in line:
    #LINE IS AN INDIVIDUAL LINE IN THE FILE
      #if "ATOM" in line[0:6]:
    # READS IN INFORMATION IF IT IS AN ATOM AND STORES IN ATOM ARRAY
        name = line[0:6]
        idnum = line[6:11]
        char = line[11:16]
        tp = line[16:20]
        att1 = line[20:23]
        att2 =line[23:26]
        x = line[26:38]
        y = line[38:46]
        z = line[46:54]
        av = line[55:61]
        bv = line[60:66]
        n = line[67:79]
        ATOM.append(Atom(name,idnum,char,tp,att1,att2,x,y,z,av,bv,n,line))

    return


def printCentroid():
    #print("in centroid")
    #print(len(ATOM))
    x = 0.0
    y = 0.0
    z = 0.0
    maxx = 0.0
    maxy = 0.0
    maxz = 0.0
    minx = 1000.0
    miny = 1000.0
    minz = 1000.0
    count = 0
    for atom in ATOM:
        x+=atom.point.x
        y+=atom.point.y
        z+=atom.point.z
        maxx = max(maxx, atom.point.x)
        minx = min(minx, atom.point.x)
        maxy = max(maxy, atom.point.y)
        miny = min(miny, atom.point.y)
        maxz = max(maxz, atom.point.z)
        minz = min(minz, atom.point.z)
        count+=1
    
    
    x/=float(count)
    y/=float(count)
    z/=float(count)
    
    print("Centroid is:")
    
    print(x, y, z)
    
    print("Ranges are:")
    print("x:  ", minx, maxx)
    print("y:  ", miny, maxy)
    print("z:  ", minz, maxz)
    
    rangex = maxx-minx
    rangey = maxy-miny
    rangez = maxz-minz
    
    rangex*=5
    rangey*=5
    rangez*=5

    rangex=round(rangex)
    rangey=round(rangey)
    rangez=round(rangez)
    
    if (args.Even):
        if (rangex%2):
            rangex+=1
        if (rangey%2):
            rangey+=1
        if (rangez%2):
            rangez+=1
    
    rangex = int(rangex)
    rangey = int(rangey)
    rangez = int(rangez)

    print("Recommended gist input command given input file (presumed ligand or binding cavity .pdb file)")
    print("gist gridspacn 0.5 gridcntr {0:6.3f} {1:6.3f} {2:6.3f} griddim {3} {4} {5}".format(x,y,z,rangex,rangey,rangez) )
    
    return



f = open(args.structure, "r")                      
#GETS ALL THE LINES AND SAVES ALL THE LINES INTO LINES
lines = f.readlines()
#LINE COUNT 
for line in lines:
    #print(line)
    parseData(line)
       
f.close()
    
printCentroid()
    
    
