#! /usr/bin/python -t 
# _*_ coding: iso-8859-1 _*_
# Last edited on 2013-02-20 00:40:37 by stolfilocal

PROG_NAME = __file__[len(__file__)-3]
PROG_DESC = "Outputs an SVG 'roadkill' diagram for said molecule"
PROG_VERS = "1.0"

PROG_COPYRIGHT = "Copyright © 2009-07-06 by the State University of Campinas (UNICAMP)"

import sys
import os
import re
import math; from math import sqrt,sin,cos,pi
sys.path[1:0] = [ '.', '..', sys.path[0] + '/../lib', os.path.expandvars('${STOLFIHOME}/lib') ]
# sys.stderr.write(re.sub('[,]', ',\n', "%s: path = %r\n" % (PROG_NAME, sys.path)));
import rn
import mformula; from mformula import MFormula
import mformula_oxocarbon
import mformula_svg; from mformula_svg import MFormula_SVG
import mformula_args;
import argparser

PROG_HELP = \
  PROG_NAME + " \\\n" \
  + argparser.help_info_HELP + " \\\n" \
  "    > {FIGURE}.svg"

PROG_INFO = \
  "NAME\n" \
  "  " + PROG_NAME + " - " + PROG_DESC + ".\n" \
  "\n" \
  "DESCRIPTION\n" \
  "  " + PROG_HELP + ".\n" \
  "\n" \
  "DOCUMENTATION OPTIONS\n" \
  + argparser.help_info_INFO + "\n" \
  "\n" \
  "AUTHOR\n" \
  "  Created 2009-07-12 by Jorge Stolfi, IC-UNICAMP.\n" \
  "\n" \
  "MODIFICATION HISTORY\n" \
  "  2009-07-12 by J. Stolfi, IC-UNICAMP: created.\n" \
  "\n" \
  "WARRANTY\n" \
  "  " + argparser.help_info_NO_WARRANTY + "\n" \
  "\n" \
  "RIGHTS\n" \
  "  " + PROG_COPYRIGHT + ".\n" \
  "\n" \
  "  " + argparser.help_info_STANDARD_RIGHTS

style = mformula_args.parse(None); 
  

def build_formula() :
  ft = build_formula_third(); # The left third of the molecule.
  
  fm = MFormula(style);
  k0 = fm.add_subformula(ft, 0,  30, [0,0]);
  k1 = fm.add_subformula(ft, 0, 150, [0,0]);
  k2 = fm.add_subformula(ft, 0, 270, [0,0]);
  fm.add_bond(k0+4,k1+0,1.5);
  fm.add_bond(k1+4,k2+0,1.5);
  fm.add_bond(k2+4,k0+0,1.5);
  return fm;

def build_formula_third() :
  "Builds the left 1/3 of the molecule, assumed centered at [0,0]."
  fs = build_formula_sixth();
  
  ft = MFormula(style);
  k0 = ft.add_subformula(fs, 0,   0, [0,0]);
  k1 = ft.add_subformula(fs, 1, 180, [0,0]);
  ft.add_bond(k0+0,k1+0,1.5);
  ft.add_bond(k0+2,k1+2,1.0);
  return ft;

def build_formula_sixth() :
  "Builds the top half of the left third of the molecule.\n" \
  "\n" \
  "  Assumes that the molecule is centered at [0,0] and that the benzene ring" \
  " has radius R.  Tilts the oxalate carboxyls so that all single bonds have length 1."
  fs = MFormula(style);
  
  cb = fs.bond_length(1.5); # Length of central cycle (benzene) bonds.
  db = fs.bond_length(2.0); # Relative length of double bonds.
  
  # Assume that the benzene ring stays hexagonal:
  R = cb;                    # Radius of benzene ring.
  c30 = cos(30*pi/180);
  s30 = sin(30*pi/180);
  
  # Allow for slight bending of the benzene-O bonds:
  a0 = 27*pi/180; # Lift angle of Benzene-O bonds.
  sys.stderr.write("a0 = %.2f\n" % (a0*180/pi));
  ca0 = cos(a0);
  sa0 = sin(a0);
  
  # Adjust the angle of C-O in oxalates so that single bonds have length 1:
  a1 = math.asin(R*s30 + 1*sa0 - 0.5); # Lift angle of C-O bond in oxalate.
  sys.stderr.write("a1 = %.2f\n" % (a1*180/pi));
  ca1 = cos(a1);
  sa1 = sin(a1);

  # Choose direction of C=O in carbonyls to bisect the corner angle:
  a2 = (pi/2 - a1)/2.0;             # Lift angle of C=O bond in oxalate.
  sys.stderr.write("a2 = %.2f\n" % (a2*180/pi));
  ca2 = cos(a2);
  sa2 = sin(a2);
  
  p0 = rn.scale(R, [-c30, +s30]);              # Upper left carbon of benzene. 
  p1 = rn.add(p0, [-ca0, +sa0]);               # Upper left ester oxygen.
  p2 = rn.add(p1, [-ca1, -sa1]);               # Upper left oxalate carbon.
  p3 = rn.add(p2, rn.scale(db, [-ca2, +sa2])); # Upper left carboxyl oxygen. 

  sys.stderr.write("p0--p1 = %7.5f\n" % rn.dist(p0,p1));
  sys.stderr.write("p1--p2 = %7.5f\n" % rn.dist(p1,p2));
  sys.stderr.write("p2--p2' = %7.5f\n" % (2*p2[1]));
  
  k0 = fs.add_atom("C", p0, 0,0);
  k1 = fs.add_atom("O", p1, 0,0);
  k2 = fs.add_atom("C", p2, 0,0);
  k3 = fs.add_atom("O", p3, 0,0);
  fs.add_bond(k0,k1,1.0);
  fs.add_bond(k1,k2,1.0);
  fs.add_bond(k2,k3,2.0);
  return fs;
    
fm = build_formula();

svg = MFormula_SVG(fm,0);
svg.standard_elem_style_1();
svg.output_figure();