cosmo_to_s_profile_ver_1.1.1.py

import numpy as np
import csv
import os
import datetime
from decimal import Decimal
starttime = datetime.datetime.now()
num_line = int(input("how many lines needed:"))
print(type(num_line))
# -----------------------get all the name list of txt file which is not transformed to csv-------------
# Get name function
def GetFileName(dir, x):
    listName = []
    listdir = os.listdir(dir)
    for fileName in os.listdir(dir):
        if os.path.splitext(fileName)[1] == x:
            fileName = os.path.splitext(fileName)[0]
            listName.append(fileName)
    return listName


def mkdir(path):
    isExists = os.path.exists(path)

    if not isExists:
        print(path + ' is created')
        os.makedirs(path)
        return True
    else:
        # 如果目录存在则不创建，并提示目录已存在
        print(path + ' have existed')
        return False


# ---- change the cosmo file to txt file under current dictionary----#
files = os.listdir('.')
for filename in files:
    portion = os.path.splitext(filename)
    if portion[1] == ".cosmo":
        newname = portion[0] + ".txt"
        os.rename(filename, newname)

cosmo_route = os.path.abspath(os.curdir) + "\\"
csv_route_area = cosmo_route + "s-profile_area\\"
csv_route_poss = cosmo_route + "s-profile\\"
mkdir(csv_route_poss)
mkdir(csv_route_area)
# Get the name of CSV and txt
name_list_csv = GetFileName(csv_route_poss, ".csv")
name_list_txt = GetFileName(cosmo_route, ".txt")

# save the difference of CSV between txt,and the difference of TXT bewteen CSV to name_list
name_list = [i for i in name_list_csv if i not in name_list_txt] + [i for i in name_list_txt if i not in name_list_csv]
print(str(len(name_list))+" cosmo files need to be transformed to s-profile")

# ------------------------get all the name list of txt file which is not transformed to csv-------------


# open the txt for all the name in name_list
for index in range(0, len(name_list)):
    name = name_list[index]
    txt = open(cosmo_route + str(name) + ".txt", "r", encoding="utf-8")
    # create and reset the transformed space for each txt.
    list = []
    coun = []
    nseg = []
    natom = []
    x = []
    y = []
    z = []
    charg = []
    area = []
    sigm0 = []
    rn = []
    sig = []
    aeff = 7.5
    reff = (aeff / 3.141592) ** 0.5
    # read line in text line by line, create a counter for line only give # and reset its value to 0
    n = 0
    line_list = txt.readlines()
    if line_list[0].split()[0] == "Gaussian":
        print("The cosmo file is based on G09")
        area_index = 2
    else:
        print("The cosmo file is based on Tmol")
        area_index = 1
    for line in line_list:
        # find the line that contain "area"
        if line.find("area") != -1:
            if line[0] == " ":
                # only the line start with " " have the area value at 2nd position
                toarea = float(line.split()[area_index]) * 0.28002851746
        if line.find("volume") != -1:
            if line[0] == " ":
                volume = float(line.split()[1])*0.14818453429
        # if counter for # is bigger than 3, that means cosmo matrix start
        if n == 3:
            list.append(line.split())
        # if # is found once,  counter was plused one.
        if len(line) == 2:
            if line[0] == "#":
                n = n + 1
    # save all cosmo data to the list. list[x,0] = (x+1)th line nseg, list [x,1]= (x+1)th line natom,
    # list[x,2]= (x+1)th line x axis,list[x,3]= (x+1)th line y axis, list[x,4] = (x+1)th line z axis,
    #  list[x,5] = (x+1)th line charge, list[x,6] =(x+1)th line area,list [x,7] = (x+1)th line sigm0
    for data in list:
        nseg.append(int(data[0]))
        natom.append(float(data[1]))
        x.append(float(data[2]))
        y.append(float(data[3]))
        z.append(float(data[4]))
        charg.append(float(data[5]))
        area.append(float(data[6]))
        sigm0.append(float(data[7]))
        rn.append((float(data[6]) / 3.141592) ** 0.5)
    # ---------------------------------------copy from process by Prof. shimoyama-------------------------------------------#
    # creat a space to save dmn2 list long =nseg, tall = nseg
    dmn2 = [([0] * nseg[-1]) for ns in range(0, nseg[-1])]
    susi = np.zeros(shape =(nseg[-1],nseg[-1]), dtype=np.float64)
    # write dmn2 with new data.
    for jj in range(0, nseg[-1]):
        for ii in range(0, nseg[-1]):
            dmn2[jj][ii] = (x[jj] - x[ii]) ** 2 + (y[jj] - y[ii]) ** 2 + (z[jj] - z[ii]) ** 2
            aa = Decimal(dmn2[jj][ii] / (rn[ii] ** 2 + reff ** 2))
            susi[jj][ii] = (1 / np.exp(aa))


    sigm = [0 for long in range(0, nseg[-1])]
    sigm1 = [0 for long in range(0, nseg[-1])]
    sigm2 = [0 for long in range(0, nseg[-1])]
    avearea = [0 for long in range(0, nseg[-1])]
    avearea1 = [0 for long in range(0, nseg[-1])]
    avearea2 = [0 for long in range(0, nseg[-1])]
    for i in range(0, nseg[-1]):
        for j in range(0, nseg[-1]):
            sigm1[i] = sigm1[i] + (sigm0[j] * (rn[j] ** 2 * reff ** 2) / (rn[j] ** 2 + reff ** 2)) * susi[i][j]
            sigm2[i] = sigm2[i] + ((rn[j] ** 2 * reff ** 2) / (rn[j] ** 2 + reff ** 2)) * susi[i][j]
            avearea1[i] = avearea1[i] + (area[j] * (rn[j] ** 2 * reff ** 2) / (rn[j] ** 2 + reff ** 2)) * susi[i][j]
            avearea2[i] = avearea2[i] + ((rn[j] ** 2 * reff ** 2) / (rn[j] ** 2 + reff ** 2)) * susi[i][j]
        sigm[i] = sigm1[i] / sigm2[i]
        avearea[i] = avearea1[i] / avearea2[i]
    # ---------------------------------------copy from process by Prof. shimoyama-------------------------------------------#
    # save data as 25 line from -0.03 to 0.03
    numd = num_line

    # coun = [0] * 61
    coun = [0] * 51
    for p in range(0, numd):
        # sig.append((-0.03 + 0.06 / (numd - 1) * p))
        sig.append((-0.025 + 0.05 / (numd - 1) * p))
        for q in range(0, nseg[-1]):
            ff = abs(sigm[q] - sig[p])
            # if ff < 0.06 / (2 * (numd - 1)):
            if ff < 0.05 / (2 * (numd - 1)):
                coun[p] = coun[p] + area[q]
            else:
                coun[p] = coun[p]

    path = csv_route_area + str(name) + ".csv"
    path_x = csv_route_poss + str(name) + ".csv"
    with open(path, "w", newline="") as f:
        csv_write = csv.writer(f)
        for g in range(0, numd):
            data_row = [sig[g], coun[g]]
            csv_write.writerow(data_row)
        csv_write.writerow(["area", toarea])
        csv_write.writerow(["volume",volume])
    with open(path_x, "w", newline="") as f:
        csv_write = csv.writer(f)
        for g in range(0, numd):
            data_row = [sig[g], (coun[g] / toarea)]
            csv_write.writerow(data_row)
        csv_write.writerow(["area", toarea])
        csv_write.writerow(["volume", volume])
    print(str(index+1) + "/" + str(len(name_list))+ " cosmo files are transformed to s-profile")
    print(name_list[index] + " is done\n")
endtime= datetime.datetime.now()
time_use = endtime - starttime
print("all " + str(len(name_list)) + " cosmo files have been transformed to s-profile")
print("totaly "+str(time_use) + " are used")
os.system("pause")