#import functools
import gui
import os
import shutil
import time
import olx
import olex
import olexex
import olex_core
import sys
import subprocess
import signal
import time
from olexFunctions import OV
#import OlexVFS
#from PIL import ImageDraw, Image
from ImageTools import IT
#from PilTools import timage
from cctbx import adptbx
from cctbx.array_family import flex
from decors import run_with_bitmap
import numpy as np

ELEMENTS = {
  1: "H", 2: "He", 3: "Li", 4: "Be", 5: "B", 6: "C", 7: "N", 8: "O", 9: "F", 10: "Ne",
  11: "Na", 12: "Mg", 13: "Al", 14: "Si", 15: "P", 16: "S", 17: "Cl", 18: "Ar",
  19: "K", 20: "Ca", 21: "Sc", 22: "Ti", 23: "V", 24: "Cr", 25: "Mn", 26: "Fe",
  27: "Co", 28: "Ni", 29: "Cu", 30: "Zn", 31: "Ga", 32: "Ge", 33: "As", 34: "Se",
  35: "Br", 36: "Kr", 37: "Rb", 38: "Sr", 39: "Y", 40: "Zr", 41: "Nb", 42: "Mo",
  43: "Tc", 44: "Ru", 45: "Rh", 46: "Pd", 47: "Ag", 48: "Cd", 49: "In", 50: "Sn",
  51: "Sb", 52: "Te", 53: "I", 54: "Xe", 55: "Cs", 56: "Ba", 57: "La", 58: "Ce",
  59: "Pr", 60: "Nd", 61: "Pm", 62: "Sm", 63: "Eu", 64: "Gd", 65: "Tb", 66: "Dy",
  67: "Ho", 68: "Er", 69: "Tm", 70: "Yb", 71: "Lu", 72: "Hf", 73: "Ta", 74: "W",
  75: "Re", 76: "Os", 77: "Ir", 78: "Pt", 79: "Au", 80: "Hg", 81: "Tl", 82: "Pb",
  83: "Bi", 84: "Po", 85: "At", 86: "Rn", 87: "Fr", 88: "Ra", 89: "Ac", 90: "Th",
  91: "Pa", 92: "U", 93: "Np", 94: "Pu", 95: "Am", 96: "Cm", 97: "Bk", 98: "Cf",
  99: "Es", 100: "Fm", 101: "Md", 102: "No", 103: "Lr", 104: "Rf", 105: "Db",
  106: "Sg", 107: "Bh", 108: "Hs", 109: "Mt", 110: "Ds", 111: "Rg", 112: "Cn",
  113: "Nh", 114: "Fl", 115: "Mc", 116: "Lv", 117: "Ts", 118: "Og"
}

ELEMENTS_BY_SYMBOL = {symbol: z for z, symbol in ELEMENTS.items()}

_interrupt_requested = False
_interrupt_handlers_installed = False
_interrupt_ignore_until = 0.0

def _signal_interrupt_handler(signum, frame):
  global _interrupt_requested
  _interrupt_requested = True

def _install_interrupt_handlers_once():
  global _interrupt_handlers_installed
  if _interrupt_handlers_installed:
    return
  for sig_name in ("SIGINT", "SIGTERM"):
    sig = getattr(signal, sig_name, None)
    if sig is None:
      continue
    try:
      signal.signal(sig, _signal_interrupt_handler)
    except Exception:
      pass
  _interrupt_handlers_installed = True

def reset_interrupt_request():
  global _interrupt_requested
  _interrupt_requested = False

def consume_interrupt_request(cooldown_s=0.6):
  global _interrupt_ignore_until
  reset_interrupt_request()
  _interrupt_ignore_until = time.time() + max(0.0, float(cooldown_s))
  if sys.platform[:3] == 'win':
    try:
      import ctypes
      ctypes.windll.user32.GetAsyncKeyState(0x13)
    except Exception:
      pass

def is_pause_break_pressed():
  _install_interrupt_handlers_once()
  if time.time() < _interrupt_ignore_until:
    return False
  if sys.platform[:3] != 'win':
    return _interrupt_requested
  try:
    import ctypes
    return bool(ctypes.windll.user32.GetAsyncKeyState(0x13) & 0x1)
  except Exception:
    return False

def terminate_process_tree(process):
  if process is None or process.poll() is not None:
    return
  try:
    if sys.platform[:3] == 'win':
      subprocess.run(
        ["taskkill", "/PID", str(process.pid), "/T", "/F"],
        stdout=subprocess.DEVNULL,
        stderr=subprocess.DEVNULL,
        check=False,
      )
    else:
      try:
        os.killpg(os.getpgid(process.pid), signal.SIGTERM)
      except Exception:
        process.terminate()
  except Exception:
    try:
      process.kill()
    except Exception:
      pass

def _get_occ_fchk_candidates_from_toml(toml_path):
  base_dir = os.path.dirname(toml_path)
  base_name = os.path.splitext(os.path.basename(toml_path))[0]
  return [
    os.path.join(base_dir, base_name + ".owf.fchk"),
    os.path.join(base_dir, "experimental.owf.fchk"),
  ]

def _collect_occ_fchk_outputs(wfn_file, folder, name):
  if isinstance(wfn_file, list):
    found = []
    for toml_path in wfn_file:
      if not str(toml_path).lower().endswith(".toml"):
        continue
      candidates = _get_occ_fchk_candidates_from_toml(toml_path)
      for candidate in candidates:
        if os.path.exists(candidate):
          found.append(candidate)
          break
    return found

  wfn_job_dir = os.path.join(folder, "olex2", "Wfn_job")
  fallback_candidates = [
    os.path.join(wfn_job_dir, name + ".owf.fchk"),
    os.path.join(wfn_job_dir, "experimental.owf.fchk"),
  ]
  if str(wfn_file).lower().endswith(".toml"):
    fallback_candidates = _get_occ_fchk_candidates_from_toml(wfn_file) + fallback_candidates

  for candidate in fallback_candidates:
    if os.path.exists(candidate):
      return [candidate]
  return []

def _copy_occ_fchk_outputs_to_workspace(fchk_files, folder, name):
  for index, source in enumerate(fchk_files, start=1):
    source_dir = os.path.basename(os.path.dirname(source))
    if source_dir.startswith("Part_"):
      suffix = "_part_" + source_dir.split("Part_", 1)[1]
    elif len(fchk_files) > 1:
      suffix = "_part_%d" % index
    else:
      suffix = ""
    destination = os.path.join(folder, name + suffix + ".owf.fchk")
    shutil.copy(source, destination)

def scrub(cmd):
  log = gui.tools.LogListen()
  olex.m(cmd)
  return log.endListen()

def write_merged_hkl():
  folder = OV.FilePath()
  name = OV.ModelSrc()
  fn = os.path.join(folder, name + "_merged.hkl")
  print("Writing " + fn)
  from cctbx_olex_adapter import OlexCctbxAdapter
  cctbx_adaptor = OlexCctbxAdapter()
  with open(fn, "w") as out:
    f_sq_obs = cctbx_adaptor.reflections.f_sq_obs_merged
    f_sq_obs.export_as_shelx_hklf(out, normalise_if_format_overflow=True)
  print("Done!")
OV.registerFunction(write_merged_hkl, False, "NoSpherA2")

@run_with_bitmap('Partitioning')
def cuqct_tsc(wfn_file, cif, groups, hkl_file=None, save_k_pts=False, read_k_pts=False):
  basis_name = OV.GetParam('snum.NoSpherA2.basis_name')
  folder = OV.FilePath()
  name = OV.ModelSrc()
  final_log_name = name + ".partitionlog"
  if not isinstance(wfn_file, list):
    gui.get_default_notification(
        txt="Calculating .tsc file from Wavefunction <b>%s</b>..."%os.path.basename(wfn_file),
        txt_col='black_text')
    if OV.HasGUI():
      OV.UpdateHtml()
      olx.Refresh()
  ncpus = OV.GetParam('snum.NoSpherA2.ncpus')
  if os.path.isfile(os.path.join(folder, final_log_name)):
    shutil.move(os.path.join(folder, final_log_name), os.path.join(folder, final_log_name + "_old"))
  args = []
  NoSpherA2 = OV.GetVar("NoSpherA2")
  args.append(NoSpherA2)
  if software() == "SALTED":
    salted_model_dir = OV.GetParam('snum.NoSpherA2.selected_salted_model')
    args.append("-SALTED")
    args.append(salted_model_dir)
  if hkl_file is not None:
    if not os.path.exists(hkl_file):
      from cctbx_olex_adapter import OlexCctbxAdapter
      cctbx_adaptor = OlexCctbxAdapter()
      f_sq_obs = cctbx_adaptor.reflections.f_sq_obs_merged
      with open(hkl_file, "w") as out:
        f_sq_obs = f_sq_obs.complete_array(d_min_tolerance=0.01, d_min=f_sq_obs.d_min() * 0.95, d_max=f_sq_obs.d_max_min()[0], new_data_value=-1, new_sigmas_value=-1)
        f_sq_obs.export_as_shelx_hklf(out, normalise_if_format_overflow=True)
  if (int(ncpus) >= 1):
    args.append('-cpus')
    args.append(ncpus)
  if OV.GetParam('snum.NoSpherA2.NoSpherA2_debug'):
    args.append('-v')
  if OV.GetParam('snum.NoSpherA2.NoSpherA2_ED'):
    args.append('-ED')
  else:
    wavelength = float(olx.xf.exptl.Radiation())
    if wavelength < 0.1:
      args.append("-ED")
  if (OV.GetParam('snum.NoSpherA2.becke_accuracy') != "Normal"):
    args.append('-acc')
    if (OV.GetParam('snum.NoSpherA2.becke_accuracy') == "Low"):
      args.append('1')
    elif (OV.GetParam('snum.NoSpherA2.becke_accuracy') == "High"):
      args.append('3')
    elif (OV.GetParam('snum.NoSpherA2.becke_accuracy') == "Max"):
      args.append('4')
  if(save_k_pts):
    args.append('-skpts')
  if(read_k_pts):
    args.append('-rkpts')
  soft = software()
  if "xTB" in soft or "pTB" in soft:
    args.append("-ECP")
    mode = OV.GetParam('snum.NoSpherA2.NoSpherA2_ECP')
    if "xTB" in soft:
      args.append(str(2))
    elif "pTB" in soft:
      args.append(str(3))
  elif "ECP" in basis_name:
    args.append("-ECP")
    mode = OV.GetParam('snum.NoSpherA2.NoSpherA2_ECP')
    args.append(str(mode))
  if OV.GetParam('snum.NoSpherA2.NoSpherA2_Partition') == "RI-Fit":
    if is_orca_new():
      auxiliary_basis = OV.GetParam('snum.NoSpherA2.auxiliary_basis')
      args.append("-RI_FIT")
      args.append(auxiliary_basis)
      if auxiliary_basis == "auto" or auxiliary_basis == "auto_aux":
        args.append(OV.GetParam('snum.NoSpherA2.auxiliary_basis_beta'))
      mem = float(OV.GetParam('snum.NoSpherA2.mem')) * 1000 #Mem in MB
      args.append("-mem")
      args.append(str(mem))
    else:
      print("WARNING! RI-FIT currently only works with ORCA newer than Version 5.0! Please update your ORCA version!")
      print("Performing refinement without RI fit!")
  elif OV.GetParam('snum.NoSpherA2.NoSpherA2_Partition') == "TFVC":
    args.append("-TFVC")
  elif OV.GetParam('snum.NoSpherA2.NoSpherA2_Partition') == "Becke":
    args.append("-Becke")
  elif OV.GetParam('snum.NoSpherA2.NoSpherA2_Partition') == "MBIS":
    args.append("-MBIS")
  elif OV.GetParam('snum.NoSpherA2.NoSpherA2_Partition') == "EMBIS":
    args.append("-EMBIS")

  olex_refinement_model = OV.GetRefinementModel(False)

  if olex_refinement_model['hklf']['value'] >= 5:
    try:
      src = OV.HKLSrc()
      cmd = "HKLF5 -e '%s'" %src
      res = scrub(cmd)
      if "HKLF5 file is expected" in " ".join(res):
        pass
      elif "negative batch numbers" in " ".join(res):
        pass
      else:
        args.append("-twin")
        for i in res[4].split() + res[6].split() + res[8].split():
          args.append(str(float(i)))
    except:
      print("There is a problem with the HKLF5 file...")
  elif 'twin' in olex_refinement_model:
    c = olex_refinement_model['twin']['matrix']
    args.append("-twin")
    for row in c:
      for el in row:
        args.append(str(float(el)))
  d2 = olex_core.GetHklStat()
  args.append("-hkl_min_max")
  args.append(str(d2['FileMinIndices'][0]))
  args.append(str(d2['FileMaxIndices'][0]))
  args.append(str(d2['FileMinIndices'][1]))
  args.append(str(d2['FileMaxIndices'][1]))
  args.append(str(d2['FileMinIndices'][2]))
  args.append(str(d2['FileMaxIndices'][2]))
  #shel = olx.Ins('SHEL')
  #omit = olx.Ins('OMIT')
  #d_min = d2['MinD']
  #if shel != "n/a":
  #  d = float(shel.split()[-1])
  #  if d > d_min:
  #    d_min = d
  #if omit != "n/a":
  #  from cctbx import uctbx
  #  d = uctbx.two_theta_as_d(float(omit.split()[-1]),float(olx.xf.exptl.Radiation()),deg=True)
  #  if d > d_min:
  #    d_min = d
  #args.append("-dmin")
  #args.append(str(d_min * 0.95))
  if OV.IsEDRefinement():
    try:
      top_up_d = float(OV.GetACI().EDI.get_stored_param("refinement.top_up_d"))
      if top_up_d < d2['MinD']:
        args.append("-dmin")
        args.append(str(top_up_d))
    except:
      print("Failed to set dmin for TSC")
  if isinstance(wfn_file, list):
    if isinstance(cif, list):
      args.append("-cmtc")
      if len(wfn_file) != len(groups) or len(wfn_file) != len(groups):
        print("Inconsistant size of parameters! ERROR!")
        return
      for i in range(len(wfn_file)):
        args.append(wfn_file[i])
        args.append(cif[i])
        for j in range(len(groups[i])):
          groups[i][j] = str(groups[i][j])
        args.append(','.join(groups[i]))
    else:
      args.append("-cif")
      args.append(cif)
      args.append("-mtc")
      if len(wfn_file) != len(groups):
        print("Inconsistant size of parameters! ERROR!")
        return
      for i in range(len(wfn_file)):
        args.append(wfn_file[i])
        for j in range(len(groups[i])):
          groups[i][j] = str(groups[i][j])
        args.append(','.join(groups[i]))
    if soft == "Hybrid":
      args.append("-mtc_mult")
      for i in range(1, min(6, len(groups) + 1)):
        m =  OV.GetParam('snum.NoSpherA2.Hybrid.multiplicity_Part%d' % i)
        if m is None or m == 'None':
          m = 1
        args.append(str(m))
      args.append("-mtc_charge")
      for i in range(1, min(6, len(groups) + 1)):
        c = int(OV.GetParam('snum.NoSpherA2.Hybrid.charge_Part%d' % i))
        if c < 0:
          args.append("n"+str(c))
        else:
          args.append(str(c))
      args.append("-mtc_ECP")
      for i in range(1, min(6, len(groups) + 1)):
        ECP_m_C = 0
        sftw = OV.GetParam('snum.NoSpherA2.Hybrid.software_Part%d' % i).lstrip()
        if sftw == "xTB":
          ECP_m_C = 2
        elif sftw == "pTB":
          ECP_m_C = 3
        args.append(str(ECP_m_C))
    else:
      args.append("-mtc_mult")
      for i in range(len(groups)):
        m =  OV.GetParam('snum.NoSpherA2.muliplicity')
        if m is None or m == 'None':
          m = 1
        args.append(str(m))
      args.append("-mtc_charge")
      for i in range(len(groups)):
        c = int(OV.GetParam('snum.NoSpherA2.charge'))
        if c < 0:
          args.append("n"+str(c))
        else:
          args.append(str(c))
      args.append("-mtc_ECP")
      for i in range(len(groups)):
        if soft == "xTB":
          args.append("2")
        elif soft == "pTB":
          args.append("3")
        elif "ECP" in basis_name:
          args.append("1")
        else:
          args.append("0")

    if any(".xyz" in f for f in wfn_file):
      Cations = OV.GetParam('snum.NoSpherA2.Thakkar_Cations')
      if Cations != "" and Cations is not None:
        args.append("-Cations")
        args.append(Cations)
      Anions = OV.GetParam('snum.NoSpherA2.Thakkar_Anions')
      if Anions != "" and Anions is not None:
        args.append("-Anions")
        args.append(Anions)
  else:
    if soft == "OCC":
      args.append("-occ")
      args.append(wfn_file)
    else:
      args.append("-wfn")
      args.append(wfn_file)
    args.append("-cif")
    args.append(cif)
    args.append("-mult")
    m = OV.GetParam('snum.NoSpherA2.multiplicity')
    if m == 0:
      m = 1
    if m is None:
      m = 1
    args.append(str(m))
    args.append("-charge")
    c = OV.GetParam('snum.NoSpherA2.charge')
    args.append(str(c))
    if(groups[0] != -1000):
      args.append('-group')
      for i in range(len(groups)):
        args.append(groups[i])
    if ".xyz" in wfn_file:
      Cations = OV.GetParam('snum.NoSpherA2.Thakkar_Cations')
      if Cations != "" and Cations is not None:
        args.append("-Cations")
        args.append(Cations)
      Anions = OV.GetParam('snum.NoSpherA2.Thakkar_Anions')
      if Anions != "" and Anions is not None:
        args.append("-Anions")
        args.append(Anions)

  env = None
  if soft == "OCC":
    env = dict(os.environ)
    env["OCC_DATA_PATH"] = os.path.join(os.path.dirname(NoSpherA2), "occ", "share")

  p = None
  if sys.platform[:3] == 'win':
    startinfo = subprocess.STARTUPINFO()
    startinfo.dwFlags |= subprocess.STARTF_USESHOWWINDOW
    startinfo.wShowWindow = 0
    creationflags = getattr(subprocess, "CREATE_NO_WINDOW", 0)
    p = subprocess.Popen(args, startupinfo=startinfo, creationflags=creationflags, env=env)
  else:
    p = subprocess.Popen(args, env=env)

  out_fn = "NoSpherA2.log"
  occ_log_fn = "NoSpherA2_OCC.log"
  wfn_log_fn = name + ".wfnlog"

  tries = 0
  while not os.path.exists(out_fn):
    if OV.GetVar("stop_current_process") or is_pause_break_pressed():
      print("Partitioning aborted by INTERRUPT before log creation!")
      terminate_process_tree(p)
      OV.SetVar("stop_current_process", False)
      consume_interrupt_request()
      OV.SetVar('NoSpherA2-Error', "Wfn-Interrupted")
      raise NameError('NoSpherA2/partitioning aborted by user.')
    if p.poll() is None:
      time.sleep(0.2)
      tries += 1
      if tries >= 10:
        if "python" in args[2] and tries <=10:
          continue
        print("Failed to locate the output file")
        OV.SetVar('NoSpherA2-Error',"NoSpherA2-Output not found!")
        raise NameError('NoSpherA2-Output not found!')
    else:
      print("process ended before the output file was found")
      OV.SetVar('NoSpherA2-Error',"NoSpherA2 ended unexpectedly!")
      raise NameError('NoSpherA2 ended unexpectedly!')
  with open(out_fn, "r") as stdout:
    occ_stdout = None
    if os.path.exists(occ_log_fn):
      occ_stdout = open(occ_log_fn, "r")
    try:
      aborted_by_user = False
      while p.poll() is None:
        x = stdout.read()
        if x:
          print(x, end='')
        if occ_stdout:
          y = occ_stdout.read()
          if y:
            print(y, end='')
        if OV.GetVar("stop_current_process") or is_pause_break_pressed():
          print("Calculation aborted by INTERRUPT!")
          terminate_process_tree(p)
          OV.SetVar("stop_current_process", False)
          consume_interrupt_request()
          aborted_by_user = True
          break
        else:
          time.sleep(0.1)
    finally:
      if occ_stdout:
        occ_stdout.close()
  if aborted_by_user:
    OV.SetVar('NoSpherA2-Error', "Wfn-Interrupted")
    raise NameError('NoSpherA2/partitioning aborted by user.')

  sucess = False
  shutil.move(out_fn, final_log_name)
  if os.path.exists(occ_log_fn):
    shutil.move(occ_log_fn, wfn_log_fn)
  with open(final_log_name, "r") as f:
    lines = f.readlines()
    if "Writing tsc file...  ... done!\n" in lines or "Writing tsc file...\n" in lines:
      sucess = True

  if sucess:
    if soft == "OCC":
      occ_fchk_files = _collect_occ_fchk_outputs(wfn_file, folder, name)
      if not occ_fchk_files:
        OV.SetVar('NoSpherA2-Error', "NoFchk")
        raise NameError('OCC run did not generate an fchk file during cuQCT call!')
      _copy_occ_fchk_outputs_to_workspace(occ_fchk_files, folder, name)

      tsc_candidates = [
        os.path.join(folder, name + ".tsc"),
        os.path.join(folder, name + ".tscb"),
        os.path.join(folder, "experimental.tsc"),
        os.path.join(folder, "experimental.tscb"),
      ]
      has_tsc = any(os.path.exists(f) for f in tsc_candidates)
      if not has_tsc:
        OV.SetVar('NoSpherA2-Error', "NoTsc")
        raise NameError('OCC run did not generate a .tsc/.tscb file during cuQCT call!')
    print("\n.tsc calculation ended!")
  else:
    print ("\nERROR during tsc calculation!")
    raise NameError('NoSpherA2-Output not complete!')

def get_nmo():
  line = ""
  if not os.path.isfile(os.path.join(OV.FilePath(),OV.ModelSrc()+".wfn")):
    part_log_path =  os.path.join(OV.FilePath(),OV.ModelSrc()+".partitionlog")
    if os.path.exists(part_log_path):
      with open(part_log_path) as partlog:
        while "Number of MOs:" not in line:
          line = partlog.readline()
      values = line.split(":")
      return int(values[2])
    else:
      return -1
  else:
    with open(os.path.join(OV.FilePath(),OV.ModelSrc()+".wfn")) as wfn:
      while "MOL ORBITAL" not in line:
        line = wfn.readline()
    values = line.split()
    return int(values[1])

def get_ncen():
  if not os.path.isfile(os.path.join(OV.FilePath(),OV.ModelSrc()+".wfn")):
    return -1
  wfn = open(os.path.join(OV.FilePath(),OV.ModelSrc()+".wfn"))
  line = ""
  while "MOL ORBITAL" not in line:
    line = wfn.readline()
  values = line.split()
  return values[6]

OV.registerFunction(get_nmo, False, 'NoSpherA2')
OV.registerFunction(get_ncen, False, 'NoSpherA2')

@run_with_bitmap("Combining .tsc")
def combine_tscs(match_phrase="_part_", no_check=False):
  gui.get_default_notification(txt="Combining .tsc files", txt_col='black_text')
  args = []
  NSP2_exe = OV.GetVar("NoSpherA2")
  args.append(NSP2_exe)
  if (no_check):
    args.append("-merge_nocheck")
  else:
    args.append("-merge")

  print("Combinging the .tsc files of disorder parts... Please wait...")
  sfc_name = OV.ModelSrc()

  tsc_dst = os.path.join(OV.FilePath(), sfc_name + "_total.tsc")
  if os.path.exists(tsc_dst):
    backup = os.path.join(OV.FilePath(), "tsc_backup")
    if not os.path.exists(backup):
      os.mkdir(backup)
    i = 1
    while (os.path.exists(os.path.join(backup,sfc_name + ".tsc") + "_%d"%i)):
      i = i + 1
    try:
      shutil.move(tsc_dst,os.path.join(backup,sfc_name + ".tsc") + "_%d"%i)
    except:
      pass

  if match_phrase != "":
    from os import walk
    _, _, filenames = next(walk(OV.FilePath()))
    for f in filenames:
      if match_phrase in f and (".tsc" in f or ".tscb" in f) :
        args.append(os.path.join(OV.FilePath(), f))
  else:
    print("ERROR! Please make sure threre is a match phrase to look for tscs!")
    return False
  startinfo = None

  if sys.platform[:3] == 'win':
    startinfo = subprocess.STARTUPINFO()
    startinfo.dwFlags |= subprocess.STARTF_USESHOWWINDOW
    startinfo.wShowWindow = subprocess.SW_HIDE
    creationflags = getattr(subprocess, "CREATE_NO_WINDOW", 0)

  if startinfo is None:
    with subprocess.Popen(args, stdout=subprocess.PIPE) as p:
      for c in iter(lambda: p.stdout.read(1), b''):
        string = c.decode()
        sys.stdout.write(string)
        sys.stdout.flush()
  else:
    with subprocess.Popen(args, stdout=subprocess.PIPE, startupinfo=startinfo, creationflags=creationflags) as p:
      for c in iter(lambda: p.stdout.read(1), b''):
        string = c.decode()
        sys.stdout.write(string)
        sys.stdout.flush()

  if os.path.exists("combined.tsc"):
    tsc_dst = sfc_name + "_total.tsc"
    shutil.move("combined.tsc", tsc_dst)
  elif os.path.exists("combined.tscb"):
    tsc_dst = sfc_name + "_total.tscb"
    shutil.move("combined.tscb", tsc_dst)
  elif os.path.exists("experimental.tscb"):
    tsc_dst = sfc_name + "_total.tscb"
    shutil.move("experimental.tscb", tsc_dst)

  try:
    OV.SetParam('snum.NoSpherA2.file', tsc_dst)
    OV.SetControlValue('SNUM_REFINEMENT_NSFF_TSC_FILE', os.path.basename(tsc_dst))
  except:
    pass
  if OV.HasGUI():
    OV.UpdateHtml()
  return True

OV.registerFunction(combine_tscs, False, 'NoSpherA2')

def deal_with_parts():
  parts = OV.ListParts()
  if not parts:
    return
  grouped_parts = read_disorder_groups()
  if grouped_parts == []:
    groups = []
    for part in parts:
      if part == 0:
        continue
      groups.append(["0",str(part)])
      olex.m("showp 0 %s" %part)
      fn = "%s_part_%s.xyz" %(OV.ModelSrc(), part)
      olx.File(fn,p=8)
    olex.m("showp")
    return list(parts), groups
  else:
    result = []
    groups = []
    longest = 0
    for i in range(len(grouped_parts)):
      if len(grouped_parts[i]) > longest:
        longest = i

    for i in range(len(grouped_parts[longest])):
      command = "showp 0"
      groups.append(["0"])
      result.append(i+1)
      for j in range(len(grouped_parts)):
        if i >= len(grouped_parts[j]):
          command += " %d"%grouped_parts[j][0]
          groups[i].append(str(grouped_parts[j][0]))
        else:
          command += " %d"%grouped_parts[j][i]
          groups[i].append(str(grouped_parts[j][i]))
      olex.m(command)
      fn = "%s_part_%s.xyz" %(OV.ModelSrc(), i+1)
      olx.File(fn,p=8)
    # I will return only a simple sequence which maps onto the calculations to be done
    olex.m("showp")
    return result, groups
OV.registerFunction(deal_with_parts, False, 'NoSpherA2')

def check_for_matching_wfn():
  p = OV.FilePath()
  name = OV.ModelSrc()
  wfn = os.path.join(p, name + '.wfn')
  wfx = os.path.join(p, name + '.wfx')
  gbw = os.path.join(p, name + '.gbw')
  ffn = os.path.join(p, name + '.ffn')
  xtb = os.path.join(p, name + '.xtb')
  if os.path.exists(wfn):
    return True
  elif os.path.exists(wfx):
    return True
  elif os.path.exists(gbw):
    return True
  elif os.path.exists(ffn):
    return True
  elif os.path.exists(xtb):
    return True
  else:
    return False
OV.registerFunction(check_for_matching_wfn, False, 'NoSpherA2')

def read_disorder_groups():
  inp = OV.GetParam('snum.NoSpherA2.Disorder_Groups')
  if inp is None or inp == "":
    return []
  groups = inp.split(';')
  result = []
  for i in range(len(groups)):
    result.append([])
    for part in groups[i].split(','):
      if '-' in part:
        if len(part) > 2:
          a, b = part.split('-')
          result[i].extend(list(range(int(a),int(b)+1)))
        else:
          result[i].append(int(part))
      else:
        result[i].append(int(part))
  return result
OV.registerFunction(read_disorder_groups, False, 'NoSpherA2')

def is_disordered():
  parts = OV.ListParts()
  if not parts:
    return False
  else:
    return True
OV.registerFunction(is_disordered, False, 'NoSpherA2')

def software():
  t = OV.GetParam('snum.NoSpherA2.source')
  t = t.lstrip().rstrip()
  return t

def org_min():
  OV.SetParam('snum.NoSpherA2.basis_name',"3-21G")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Low")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"SloppySCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"EasyConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.01")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.6")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',False)
  OV.SetParam('snum.NoSpherA2.full_HAR',False)
  olex.m("html.Update()")
OV.registerFunction(org_min, False, "NoSpherA2")

def org_small():
  OV.SetParam('snum.NoSpherA2.basis_name',"def2-SVP")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Low")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"NoSpherA2SCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"EasyConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.001")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.6")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',False)
  OV.SetParam('snum.NoSpherA2.full_HAR',False)
  olex.m("html.Update()")
OV.registerFunction(org_small, False, "NoSpherA2")

def org_final():
  OV.SetParam('snum.NoSpherA2.basis_name',"cc-pVTZ")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Normal")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"StrongSCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"EasyConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.0001")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.6")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',False)
  OV.SetParam('snum.NoSpherA2.full_HAR',True)
  olex.m("html.Update()")
OV.registerFunction(org_final, False, "NoSpherA2")

def light_min():
  OV.SetParam('snum.NoSpherA2.basis_name',"3-21G")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Low")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"SloppySCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"SlowConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.01")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.85")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',False)
  OV.SetParam('snum.NoSpherA2.full_HAR',False)
  olex.m("html.Update()")
OV.registerFunction(light_min, False, "NoSpherA2")

def light_small():
  OV.SetParam('snum.NoSpherA2.basis_name',"def2-SVP")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Low")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"NoSpherA2SCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"SlowConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.001")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.85")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',False)
  OV.SetParam('snum.NoSpherA2.full_HAR',False)
  olex.m("html.Update()")
OV.registerFunction(light_small, False, "NoSpherA2")

def light_final():
  OV.SetParam('snum.NoSpherA2.basis_name',"def2-TZVP")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Normal")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"StrongSCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"SlowConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.0001")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.85")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',False)
  OV.SetParam('snum.NoSpherA2.full_HAR',True)
  olex.m("html.Update()")
OV.registerFunction(light_final, False, "NoSpherA2")

def heavy_min():
  OV.SetParam('snum.NoSpherA2.basis_name',"x2c-SVP")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Low")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"SloppySCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"SlowConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.01")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.85")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',True)
  OV.SetParam('snum.NoSpherA2.full_HAR',False)
  olex.m("html.Update()")
OV.registerFunction(heavy_min, False, "NoSpherA2")

def heavy_small():
  OV.SetParam('snum.NoSpherA2.basis_name',"jorge-DZP-DKH")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Low")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"NoSpherA2SCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"SlowConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.001")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.85")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',True)
  OV.SetParam('snum.NoSpherA2.full_HAR',False)
  olex.m("html.Update()")
OV.registerFunction(heavy_small, False, "NoSpherA2")

def heavy_final():
  OV.SetParam('snum.NoSpherA2.basis_name',"x2c-TZVP")
  if "Tonto" in software():
    OV.SetParam('snum.NoSpherA2.method', "B3LYP")
  elif "ORCA" in software():
    if is_orca_new():
      OV.SetParam('snum.NoSpherA2.method', "r2SCAN")
    else:
      OV.SetParam('snum.NoSpherA2.method', "PBE")
  else:
    OV.SetParam('snum.NoSpherA2.method', "PBE")
  OV.SetParam('snum.NoSpherA2.becke_accuracy',"Normal")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Conv',"StrongSCF")
  OV.SetParam('snum.NoSpherA2.ORCA_SCF_Strategy',"SlowConv")
  OV.SetParam('snum.NoSpherA2.cluster_radius',0)
  OV.SetParam('snum.NoSpherA2.DIIS',"0.0001")
  OV.SetParam('snum.NoSpherA2.pySCF_Damping',"0.85")
  OV.SetParam('snum.NoSpherA2.ORCA_Solvation',"Water")
  OV.SetParam('snum.NoSpherA2.Relativistic',True)
  OV.SetParam('snum.NoSpherA2.full_HAR',True)
  olex.m("html.Update()")
OV.registerFunction(heavy_final, False, "NoSpherA2")

def make_quick_button_gui():
  import olexex

  max_Z, heaviest_element_in_formula = olexex.FindZOfHeaviestAtomInFormula()
  d = {}
  if max_Z < 10:
    d.setdefault('type', 'org')
    d.setdefault('description', 'Organic Molecule (Z &lt; 10) ')
  elif max_Z < 36:
    d.setdefault('type', 'light')
    d.setdefault('description', 'Light metal (Z &lt; 35) ')
  else:
    d.setdefault('type', 'heavy')
    d.setdefault('description', 'Heavy metal (Z  &gt; 35) ')

  buttons = ""
  base_col = OV.GetParam("gui.action_colour")
  for item in [("min", IT.adjust_colour(base_col, luminosity=1.0, as_format='hex'), "Test"),
               ("small", IT.adjust_colour(base_col, luminosity=0.9, as_format='hex'), "Work"),
               ("final", IT.adjust_colour(base_col, luminosity=0.8, as_format='hex'), "Final")]:
    d["qual"]=item[0]
    d["col"]=item[1]
    d["value"]=item[2]
    buttons += '''
    <td width="20%%">
      <input
        type="button"
        name="nosphera2_quick_button_%(type)s_%(qual)s"
        value="%(value)s"
        height="GetVar(HtmlComboHeight)"
        onclick="spy.NoSpherA2.%(type)s_%(qual)s()"
        bgcolor="%(col)s"
        fgcolor="#ffffff"
        fit="false"
        flat="GetVar(linkButton.flat)"
        hint=""
        disabled="false"
        custom="GetVar(custom_button)"
      >
    </td>''' % d

  d["buttons"]=buttons

  t='''
<td width="40%%" align='left'>
<b>%(description)s</b>
</td>
%(buttons)s
''' % d
  return t
OV.registerFunction(make_quick_button_gui, False, "NoSpherA2")

def calculate_number_of_electrons():
  ne = -int(OV.GetParam('snum.NoSpherA2.charge'))
  model = olexex.OlexRefinementModel()
  for sc in model._atoms:
    t = sc['type']
    if t == 'Q':
      continue
    Z = ELEMENTS_BY_SYMBOL.get(t, 0)
    ne += Z
  return ne, model

def source_is_tsc():
  source = software()
  if ".tsc" in source or ".tscb" in source:
    return True
  else:
    return False
OV.registerFunction(source_is_tsc, False, 'NoSpherA2')

def source_is_discamb():
  source = software()
  if source == OV.GetParam('user.NoSpherA2.discamb_exe'):
    return True
  else:
    return False
OV.registerFunction(source_is_discamb, False, 'NoSpherA2')

def write_precise_model_file(model = None, cov_matrix = None, annotations = None):
  from refinement import FullMatrixRefine
  from olexex import OlexRefinementModel
  table_name = ""
  if OV.IsNoSpherA2():
    table_name = str(OV.GetParam("snum.NoSpherA2.file"))
  try:
    fmr = FullMatrixRefine()
    if table_name != "":
      norm_eq = fmr.run(build_only=True, table_file_name=table_name)
    else:
      norm_eq = fmr.run(build_only=True)
    norm_eq.build_up(False)
    reparam = norm_eq.reparametrisation
    reparam.linearise()
    jac_tr = reparam.jacobian_transpose_matching_grad_fc()
    cov_matrix = flex.abs(flex.sqrt(norm_eq.covariance_matrix().matrix_packed_u_diagonal()))
    esds = jac_tr.transpose() * flex.double(cov_matrix)
    jac_tr = None
    annotations = reparam.component_annotations
  except:
    print("Could not obtain cctbx object and calculate ESDs!\n")
    return False
  f = open(OV.ModelSrc() + ".precise_model", "w")
  matrix_run = 0
  model = OlexRefinementModel()
  UC = norm_eq.xray_structure.unit_cell()
  f.write("Positions     a          b          c       and ESDs:\n")
  for atom in model._atoms:
    xyz = atom['crd'][0]
    f.write("{:5}".format(atom['label']))
    for x in range(3):
      f.write("{:12.8f}".format(xyz[x]))
    for x in range(3):
      f.write("{:12.8f}".format(esds[matrix_run + x]))
    matrix_run += 3
    has_adp = atom.get('adp')
    if has_adp is not None:
      matrix_run += 6
    else:
      matrix_run += 1
    if matrix_run < len(annotations):
      if ".C111" in annotations[matrix_run]:
        matrix_run += 25
    if matrix_run < len(annotations):
      if 'occ' in annotations[matrix_run]:
        matrix_run += 1
    if matrix_run < len(annotations):
      if 'fp' in annotations[matrix_run]:
        matrix_run += 2
    f.write("\n")
  matrix_run = 0
  f.write("\n\nADPs      (11)        (22)        (33)        (23)        (13)        (12)        Ueq     And ESDs\n")
  for atom in model._atoms:
    has_adp = atom.get('adp')
    f.write("{:5}".format(atom['label']))
    matrix_run += 3
    if has_adp is not None:
      adp = has_adp[0]
      adp = (adp[0], adp[1], adp[2], adp[5], adp[4], adp[3])
      uiso = adptbx.u_cart_as_u_iso(adp)
      adp = adptbx.u_cart_as_u_cif(UC, adp)
      for u in range(6):
        f.write("{:12.8f}".format(adp[u]))
      f.write("{:12.8f}".format(uiso))
      adp_esds = (esds[matrix_run],
                  esds[matrix_run + 1],
                  esds[matrix_run + 2],
                  esds[matrix_run + 3],
                  esds[matrix_run + 4],
                  esds[matrix_run + 5])
      u_iso_esd = adptbx.u_star_as_u_iso(UC, adp_esds)
      adp_esds = adptbx.u_star_as_u_cif(UC, adp_esds)
      for u in range(6):
        f.write("{:12.8f}".format(adp_esds[u]))
      f.write("{:12.8f}".format(u_iso_esd))
      matrix_run += 6
    else:
      for u in range(6):
        f.write("{:12s}".format(" ---"))
      f.write("{:12.8f}".format(atom['uiso'][0]))
      for u in range(6):
        f.write("{:12s}".format(" ---"))
      f.write("{:12.8f}".format(esds[matrix_run]))
      matrix_run += 1
    if matrix_run < len(annotations):
      if ".C111" in annotations[matrix_run]:
        matrix_run += 25
    if matrix_run < len(annotations):
      if 'occ' in annotations[matrix_run]:
        matrix_run += 1
    if matrix_run < len(annotations):
      if 'fdp' in annotations[matrix_run]:
        matrix_run += 2
    f.write("\n")

  connectivity_full = fmr.reparametrisation.connectivity_table
  xs = fmr.xray_structure()
  from scitbx import matrix
  import math
  from cctbx.crystal import calculate_distances
  cell_params = fmr.olx_atoms.getCell()
  cell_errors = fmr.olx_atoms.getCellErrors()
  cell_vcv = flex.pow2(matrix.diag(cell_errors).as_flex_double_matrix())
  dist_stats = {}
  dist_errs = {}
  for i in range(3):
    for j in range(i+1,3):
      if (cell_params[i] == cell_params[j] and
            cell_errors[i] == cell_errors[j] and
            cell_params[i+3] == 90 and
            cell_errors[i+3] == 0 and
            cell_params[j+3] == 90 and
            cell_errors[j+3] == 0):
        cell_vcv[i,j] = math.pow(cell_errors[i],2)
        cell_vcv[j,i] = math.pow(cell_errors[i],2)
  #Prepare the Cell Variance covariance matrix, since we need it for error propagation in distances
  cell_vcv = cell_vcv.matrix_symmetric_as_packed_u()
  sl = xs.scatterers().extract_labels()
  sf = xs.sites_frac()
  #This is VCV from refinement equations
  cm = norm_eq.covariance_matrix_and_annotations().matrix
  pm = xs.parameter_map()
  pat = connectivity_full.pair_asu_table

  # calculate the distances using the prepared information
  distances = calculate_distances(
              pat,
              sf,
              covariance_matrix=cm,
              cell_covariance_matrix=cell_vcv,
              parameter_map=pm)

  #The distances only exist once we iterate over them! Therefore build them and save them in this loop
  for i,d in enumerate(distances):
    bond = sl[d.i_seq]+"-"+sl[d.j_seq]
    dist_stats[bond] = distances.distances[i]
    var = distances.variances[i]
    if (var > 0) :
      dist_errs[bond] = math.sqrt(var)
    else:
      if abs(var) > 1E-12:
        dist_errs[bond] = math.sqrt(var)
      else:
        dist_errs[bond] = 0.0

  f.write("\n\nBondlengths and errors:\n")
  for key in dist_stats:
    f.write(f"{key} {dist_stats[key]} {dist_errs[key]}\n")

  f.close()
OV.registerFunction(write_precise_model_file, False, "NoSpherA2")

def get_tsc_file_dropdown_items():
    res = gui.GetFileListAsDropdownItems(OV.FilePath(), "tsc;tscb")
    t = res.split(";")
    if len(t[0]) == 0:
        OV.SetParam('snum.NoSpherA2.file', 'No .tsc/.tscb files found')
        return "No .tsc/.tscb files found"
    else:
        return res
OV.registerFunction(get_tsc_file_dropdown_items, False, "NoSpherA2")

@run_with_bitmap("Polarizibilities")
def calc_polarizabilities(efield = 0.005, resolution = 0.1, radius = 2.5):
  from NoSpherA2.NoSpherA2 import NoSpherA2_instance as nsp2
  from . import Wfn_Job
  pol_fol =  os.path.join("olex2", "Polarizability")
  if not os.path.exists(pol_fol):
    try:
      os.mkdir(pol_fol)
    except:
      pass
  OV.SetVar("Run_number", 0)
  wfn_job = Wfn_Job.wfn_Job(nsp2, olx.FileName(), dir="olex2\\Polarizability")
  wfn_job.write_orca_input(True, efield="00")
  wfn_job.run(copy=False)
  gbw_name = os.path.join(pol_fol, olx.FileName() + ".gbw")
  zero =  os.path.join(pol_fol, "zero.gbw")
  shutil.move(gbw_name, zero)
  wfn_job.write_orca_input(True, efield=f"x{efield}")
  wfn_job.run(copy=False)
  xp = os.path.join(pol_fol, "xp.gbw")
  shutil.move(gbw_name, xp)
  wfn_job.write_orca_input(True, efield=f"x-{efield}")
  wfn_job.run(copy=False)
  xm = os.path.join(pol_fol, "xm.gbw")
  shutil.move(gbw_name, xm)
  wfn_job.write_orca_input(True, efield=f"y{efield}")
  wfn_job.run(copy=False)
  yp = os.path.join(pol_fol, "yp.gbw")
  shutil.move(gbw_name, yp)
  wfn_job.write_orca_input(True, efield=f"y-{efield}")
  wfn_job.run(copy=False)
  ym = os.path.join(pol_fol, "ym.gbw")
  shutil.move(gbw_name, ym)
  wfn_job.write_orca_input(True, efield=f"z{efield}")
  wfn_job.run(copy=False)
  zp = os.path.join(pol_fol, "zp.gbw")
  shutil.move(gbw_name, zp)
  wfn_job.write_orca_input(True, efield=f"z-{efield}")
  wfn_job.run(copy=False)
  zm = os.path.join(pol_fol, "zm.gbw")
  shutil.move(gbw_name, zm)
  args = []
  NSP2_exe = OV.GetVar("NoSpherA2")
  args.append(NSP2_exe)
  args.append("-polarizabilities")
  args.append(zero)
  args.append(xp)
  args.append(xm)
  args.append(yp)
  args.append(ym)
  args.append(zp)
  args.append(zm)
  args.append("-e_field")
  args.append(str(efield))
  args.append("-resolution")
  args.append(str(resolution))
  args.append("-radius")
  args.append(str(radius))
  args.append("-cpus")
  args.append(OV.GetParam("snum.NoSpherA2.ncpus"))

  startinfo = None

  from subprocess import Popen, PIPE
  from sys import stdout
  if sys.platform[:3] == 'win':
    from subprocess import STARTUPINFO, STARTF_USESHOWWINDOW, SW_HIDE
    startinfo = STARTUPINFO()
    startinfo.dwFlags |= STARTF_USESHOWWINDOW
    startinfo.wShowWindow = SW_HIDE
    creationflags = getattr(subprocess, "CREATE_NO_WINDOW", 0)

  if startinfo is None:
    with Popen(args, stdout=PIPE) as p:
      for c in iter(lambda: p.stdout.read(1), b''):
        string = c.decode()
        stdout.write(string)
        stdout.flush()
  else:
    with Popen(args, stdout=PIPE, startupinfo=startinfo, creationflags=creationflags) as p:
      for c in iter(lambda: p.stdout.read(1), b''):
        string = c.decode()
        stdout.write(string)
        stdout.flush()

  NSA2_log = "NoSpherA2.log"
  if os.path.exists(NSA2_log):
    shutil.copy(NSA2_log, olx.FileName()+".polarizability")

OV.registerFunction(calc_polarizabilities, namespace="NoSpherA2")

#CHOOSE FOLDER UTILITIES! (Used in SALTED)
def setDir(phil_value) -> None:
  """Choose a directory and set the given phil value to the result.

      Args:
          None

      Returns:
          str: The chosen directory path.
      """
  try:
    a = olex.f('choosedir("Choose your data folder")')
  except:
    print("No selection made, Aborting!")
    return
  OV.SetParam(phil_value, a)
OV.registerFunction(setDir, False, "NoSpherA2")

def appendDir(phil_value) -> None:
  """Choose a directory and set teh given phil value to the result.

      Args:
          None

      Returns:
          str: The chosen directory path.
      """
  try:
    a = olex.f('choosedir("Choose your data folder")')
  except:
    print("No selection made, Aborting")
    return
  old = OV.GetParam(phil_value)
  if old is None:
    old = a
  else:
    old += ";"+a
  OV.SetParam(phil_value, old)
OV.registerFunction(appendDir, False, "NoSpherA2")

def removeDir(phil_value, item_to_remove) -> None:
  """Choose a directory and set the given phil value to the result.

      Args:
          None

      Returns:
          str: The chosen directory path.
      """
  if item_to_remove == "" or item_to_remove == "Please Select":
    return
  item_to_remove = ";" + item_to_remove
  old = OV.GetParam(phil_value)
  old = old.replace(item_to_remove,"")

  OV.SetParam(phil_value, old)
OV.registerFunction(removeDir, False, "NoSpherA2")

def toggle_full_HAR():
  if OV.GetParam('snum.NoSpherA2.full_HAR'):
    OV.SetParam('snum.NoSpherA2.full_HAR', False)
  else:
    OV.SetParam('snum.NoSpherA2.full_HAR', True)
  OV.SetItemState("h3-NoSpherA2-extras 2")
  OV.SetItemState("h3-NoSpherA2-extras 1")
OV.registerFunction(toggle_full_HAR, False, "NoSpherA2")

def cov_mat():
  with open(os.path.join(olx.FilePath(), olx.FileName()) + ".npy", "rb") as x:
    l = x.readline()
    if l != b"VCOV\n":
      print("Wrong file %s" %l)
    ann = x.readline().split()
    nf = np.load(x)
    header_string =""
    for anno in ann:
      header_string = header_string  + str(anno.decode('UTF-8'))+ ","
    sqrt_size = len(ann)
    cov_array = np.zeros((sqrt_size, sqrt_size))
    idx = 0
    for i in range(sqrt_size):
        for j in range(i, sqrt_size):
            cov_array[i, j] = nf[idx]
            cov_array[j, i] = nf[idx]
            idx += 1
    np.savetxt("numpy.csv", cov_array, delimiter=',')
    corrMat = np.corrcoef(cov_array)
    np.savetxt("./correlationMat.csv", corrMat, delimiter=",", header = header_string)

    # After calculating corrMat
    threshold = 0.75
    iu_nodiag = np.triu_indices(sqrt_size, k=1)
    strong_mask = np.abs(corrMat[iu_nodiag]) > threshold

    print(f"Correlations above {threshold}:")
    for i, j, corr in zip(iu_nodiag[0][strong_mask], iu_nodiag[1][strong_mask], corrMat[iu_nodiag][strong_mask]):
        print(f"{ann[i].decode('utf-8')} <-> {ann[j].decode('utf-8')}: {corr:.3f}")

OV.registerFunction(cov_mat, False, "NoSpherA2")

def is_orca_new(given_software = None):
  """Check if the current wavefunction software is ORCA."""
  if given_software is None:
    return software() in ["ORCA 5.0", "ORCA 6.0", "ORCA 6.1"]
  else:
    return given_software in ["ORCA 5.0", "ORCA 6.0", "ORCA 6.1"]
OV.registerFunction(is_orca_new, False, "NoSpherA2")