#!/usr/bin/env python

# eigenhist.py : converts the cor_EIG.dat file output by CA S into
#               a file of gnuplot commands to view as a histogram.

# use: python eigenhist.py cor_EIG.dat [eigenhist.dat] > gnuplotfile
#   where cor_EIG.dat is the input file from CA S,
#   eigenhist.dat, is the output doc file of percents, with keys as the eigenvalue number
#   If no output is specified, the output is named 'eigenhist.ext'
"""
cor_EIG is a text file that stores the eigenvector values. The first line
should have the number of factors on the first line: 
      1. NFAC = No. of factors 
      2. TOTAL WEIGHT = the sum of all the pixels used. Not the # of them. 
      3. TRACE = sum of the diagonal of the eigenvector array 
    There is a line for each factor. The next NFAC lines have 3 columns: 
      1.EIGENVALUE = the eigenvalue, listed largest first. 
      2.% = the % of the total interimage variance this eigenvalue is responsible for. 
      3.CUMULATIVE % = a running tally of how much interimage variance is accounted for by eigenvalues 
"""

import os,sys
from Spider.Spiderutils import writedoc

def makeGnuplotCommands(filename, P):
    nboxes = len(P)
    min = 0.2
    max = nboxes + 0.5
    gnutxt = 'set xlabel "Eigenvalue number"\n'
    gnutxt += 'set ylabel "%"\n'
    gnutxt += 'set xrange [%3.1f:%3.1f]\n' % (min, max)
    gnutxt += 'set boxwidth 0.5\n'
    gnutxt += 'plot "%s" using 1:3 title "%s" with boxes\n' % (filename, filename)
    return gnutxt

#--------------------------------------------------------------------
if __name__ == '__main__':
    length = len(sys.argv[1:])

    docfile = 'eigenhist'

    if length == 0:
        print "Usage: python eigenhist.py PREFIX_EIG_file [output_docfile] > gnuplotfile"
        sys.exit(0)
    elif length == 1:
        eigfile = sys.argv[1]
        name, ext = os.path.splitext(eigfile)
        docfile += ext
    elif length == 2:
        eigfile = sys.argv[1]
        docfile = sys.argv[2]

    # read the EIG file
    fp = open(eigfile,'r')
    B = fp.readline()   # get the number of factors from the first line
    try:
        nfact = int(B.split()[0])
    except:
        print "unable to get number of factors from %s" % eigfile
        fp.close()
        sys.exit()
    P = []
    for i in range(nfact):
        line = fp.readline()
        try:
            percent = float(line.split()[1])
            P.append(percent)
        except:
            print line
        
    fp.close()

    writedoc(docfile, columns=[P], headers=['percent'])

    print makeGnuplotCommands(docfile, P)