#!/usr/bin/env python
# -*- coding: utf-8 -*-


####################################################
# Copyright (C) 2010 Olivier Geoffroy, INOPRO
#
#	import listing file from Code Saturne and plot residuals
#
#	This program is free software: you can redistribute it and/or modify
#	it under the terms of the GNU General Public License as published by
#	the Free Software Foundation, either version 3 of the License, or
#	(at your option) any later version.
#	This program is distributed in the hope that it will be useful,
#	but WITHOUT ANY WARRANTY; without even the implied warranty of
#	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
#	GNU General Public License for more details.
#	You should have received a copy of the GNU General Public License
#	along with this program.  If not, see <http://www.gnu.org/licenses/>.
#
####################################################

import sys, os, string, pylab
from string import *


class CSlist :
	""" module for parsing listing file from code_saturne run and plot residuals"""
	def __init__(self):
		self.nbiter=[]
		self.yplusmin=[]
		self.yplusmax=[]
		self.convP=[]
		self.convVx=[]
		self.convVy=[]
		self.convVz=[]
		self.convEturb=[]
		self.convDissip=[]
		self.convphi=[]
		self.convfbarre=[]
		self.convomega=[]
		self.convR11=[]
		self.convR22=[]
		self.convR33=[]
		self.convR12=[]
		self.convR13=[]
		self.convR23=[]
		self.convTempK=[]
		self.convTempC=[]

	def load_listing(self, infile):
		""" load and parse listing """
		print "Importing data from listing "+infile
		f=open(infile,'r')
		for line in f.readlines() :
			# check for iteration
			if (line.find("TIME STEP NUMBER") > 0):
				self.nbiter.append(atoi(split(line)[5]))
			elif (line.find("Dimensionless distance   yplus") >= 0):
				self.yplusmin.append(atof(split(line)[4]))
				self.yplusmax.append(atof(split(line)[5]))
			elif (line.find("c  Pression") >= 0):
				self.convP.append(string.atof(split(line)[2]))
			elif (line.find("c  VitesseX") >= 0):
				self.convVx.append(string.atof(split(line)[2]))
			elif (line.find("c  VitesseY") >= 0):
				self.convVy.append(string.atof(split(line)[2]))
			elif (line.find("c  VitesseZ") >= 0):
				self.convVz.append(string.atof(split(line)[2]))
			elif (line.find("c  EnerTurb") >= 0):
				self.convEturb.append(string.atof(split(line)[2]))
			elif (line.find("c  Dissip") >= 0):
				self.convDissip.append(string.atof(split(line)[2]))
			elif (line.find("c  phi") >= 0):
				self.convphi.append(string.atof(split(line)[2]))
			elif (line.find("c  f_barre") >= 0):
				self.convfbarre.append(string.atof(split(line)[2]))
			elif (line.find("c  omega") >= 0):
				self.convomega.append(string.atof(split(line)[2]))
			elif (line.find("c  R11 ") >= 0):
				self.convR11.append(string.atof(split(line)[2]))
			elif (line.find("c  R22") >= 0):
				self.convR22.append(string.atof(split(line)[2]))
			elif (line.find("c  R33") >= 0):
				self.convR33.append(string.atof(split(line)[2]))
			elif (line.find("c  R12") >= 0):
				self.convR12.append(string.atof(split(line)[2]))
			elif (line.find("c  R13") >= 0):
				self.convR13.append(string.atof(split(line)[2]))
			elif (line.find("c  R23") >= 0):
				self.convR23.append(string.atof(split(line)[2]))
			elif (line.find("c  TempK") >= 0):
				self.convTempK.append(string.atof(split(line)[2]))
			elif (line.find("c  TempC") >= 0):
				self.convTempC.append(string.atof(split(line)[2]))
		f.close()
		if (len(self.nbiter) != len(self.convP)) :
			print "WARRNING: incorect dimension"
			self.nbiter.pop()


	def print_info(self):
		"""print info """
		nb=len(self.nbiter)
		print "nbr iter: ",nb," , nbr valeur convP: ",len(self.convP)," , nbr valeur convVX: ",len(self.convVx)," , nbr valeur yplus min: ",len(self.yplusmin)
		if (nb != len(self.convP)) :
			print "WARRNING: incoherence dans l'import des  données listing, le solveur à probablement divergé"
		else :
			print "import fait nbr d iteration totale = ",nb

	def export(self,fname):
		"""export data in file for post processing """
		nb=len(self.nbiter)
		print "exporting data ..."
		f=open(fname,'w')
		for i in xrange(0,nb) :
			f.write(str(self.nbiter[i])+" "+str(self.convP[i])+" "+str(self.convVx[i])+" "+str(self.convVy[i])+" "+str(self.convVz[i])+" "+str(self.yplusmin[i])+" "+str(self.yplusmax[i])+"\n")
		f.close()

	def plot_yplus(self):
		"""plot yplus values"""
		print "ploting Yplus values ..."
		pylab.plot(self.nbiter,self.yplusmin,'b-',self.nbiter,self.yplusmax,'g-')
		pylab.axhline(y=1,color='r')
		pylab.axhline(y=30,color='r')
		pylab.axhline(y=100,color='r')
		pylab.legend(("Yplus min","Yplus max"),'upper center')
		pylab.title("monitoring Yplus")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_RijSSG_TK(self):
		""" plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempK,'b+-',self.nbiter,self.convDissip,'y-',self.nbiter,self.convR11,'b--',self.nbiter,self.convR22,'g--',self.nbiter,self.convR33,'r--',self.nbiter,self.convR12,'c--',self.nbiter,self.convR13,'y--',self.nbiter,self.convR23,'k--')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","TempK","Dissip","R11","R22","R33","R12","R13","R23"),'upper center')
		pylab.title("Residus RijSSG with temp (K)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_RijSSG_TC(self):
		""" plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempC,'b+-',self.nbiter,self.convDissip,'y-',self.nbiter,self.convR11,'b--',self.nbiter,self.convR22,'g--',self.nbiter,self.convR33,'r--',self.nbiter,self.convR12,'c--',self.nbiter,self.convR13,'y--',self.nbiter,self.convR23,'k--')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","TempC","Dissip","R11","R22","R33","R12","R13","R23"),'upper center')
		pylab.title("Residus RijSSG with temp (C)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_RijSSG(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convDissip,'y-',self.nbiter,self.convR11,'b--',self.nbiter,self.convR22,'g--',self.nbiter,self.convR33,'r--',self.nbiter,self.convR12,'c--',self.nbiter,self.convR13,'y--',self.nbiter,self.convR23,'k--')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","Dissip","R11","R22","R33","R12","R13","R23"),'upper center')
		pylab.title("Residus RijSSG")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_kwSST_TK(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempK,'b+',self.nbiter,self.convEturb,'m-',self.nbiter,self.convomega,'y-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","TempK","vitesse Z","Eturb","omega"),'upper center')
		pylab.title("Residus kwSST with temp (K)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_kwSST_TC(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempC,'b+',self.nbiter,self.convEturb,'m-',self.nbiter,self.convomega,'y-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","TempC","vitesse Z","Eturb","omega"),'upper center')
		pylab.title("Residus kwSST with temp (C)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_kwSST(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convEturb,'m-',self.nbiter,self.convomega,'y-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","Eturb","omega"),'upper center')
		pylab.title("Residus kwSST")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_v2f_TK(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempK,'b+',self.nbiter,self.convEturb,'m-',self.nbiter,self.convDissip,'y-',self.nbiter,self.convphi,'k-',self.nbiter,self.convfbarre,'kx')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","TempK","Eturb","Dissip","phi","fbarre"),'upper center')
		pylab.title("Residus v2f 3D with temp (K)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_v2f_TC(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempC,'b+',self.nbiter,self.convEturb,'m-',self.nbiter,self.convDissip,'y-',self.nbiter,self.convphi,'k-',self.nbiter,self.convfbarre,'kx')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","TempC","Eturb","Dissip","phi","fbarre"),'upper center')
		pylab.title("Residus v2f with temp (C)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_v2f(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convEturb,'m-',self.nbiter,self.convDissip,'y-',self.nbiter,self.convphi,'k-',self.nbiter,self.convfbarre,'kx')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","Eturb","Dissip","phi","fbarre"),'upper center')
		pylab.title("Residus v2f")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_ke_TK(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempK,'b+',self.nbiter,self.convEturb,'m-',self.nbiter,self.convDissip,'y-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","TempK","Eturb","Dissip"),'upper center')
		pylab.title("Residus ke with temp (K)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_ke_TC(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempC,'b+',self.nbiter,self.convEturb,'m-',self.nbiter,self.convDissip,'y-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","TempC","Eturb","Dissip"),'upper center')
		pylab.title("Residus ke with temp (C)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_ke(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convEturb,'m-',self.nbiter,self.convDissip,'y-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","Eturb","Dissip"),'upper center')
		pylab.title("Residus ke")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_lam_TK(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempK,'m-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","TempK"),'upper center')
		pylab.title("Residus laminaire with temp (K)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_lam_TC(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-',self.nbiter,self.convTempC,'m-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z","TempC"),'upper center')
		pylab.title("Residus laminaire with temp (C)")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_lam(self):
		"""plot residuals"""
		pylab.semilogy(self.nbiter,self.convP,'b-',self.nbiter,self.convVx,'g-',self.nbiter,self.convVy,'r-',self.nbiter,self.convVz,'c-')
		pylab.legend(("Pression","Vitesse X","vitesse Y","vitesse Z"),'upper center')
		pylab.title("Residus laminaire")
		pylab.xlabel("iterations")
		pylab.ylabel("residus")
		pylab.grid()
		pylab.show()

	def plot_res(self):
		""" plot résiduals values"""
		if (len(self.convR11) > 0) :
			if (len(self.convTempK) > 0) :
				self.plot_RijSSG_TK()
			elif (len(self.convTempC) > 0) :
				self.plot_RijSSG_TC()
			else :
				self.plot_RijSSG()
		elif (len(self.convomega) > 0) :
			if (len(self.convTempK) > 0) :
				self.plot_kwSST_TK()
			elif (len(self.convTempC) > 0) :
				self.plot_kwSST_TC()
			else :
				self.plot_kwSST
		elif (len(self.convfbarre) > 0) :
			if (len(self.convTempK) > 0) :
				self.plot_v2f_TK()
			elif (len(self.convTempC) > 0) :
				self.plot_v2f_TC()
			else :
				self.plot_v2f()
		elif (len(self.convEturb) > 0) :
			if (len(self.convTempK) > 0) :
				self.plot_ke_TK()
			elif (len(self.convTempC) > 0) :
				self.plot_ke_TC()
			else :
				self.plot_ke()
		else :
			if (len(self.convTempK) > 0) :
				self.plot_lam_TK()
			elif (len(self.convTempC) > 0) :
				self.plot_lam_TC()
			else :
				self.plot_lam()



progname = os.path.basename(sys.argv[0])
progversion = "0.2"

if len(sys.argv) <= 1:
	print progname+" version "+progversion+" :plot residuals from Code_Saturne listing files"
	sys.exit("Usage: "+progname+" [option] listing1 listing2 ...\n   -yplus : plot yplus min and max\n   -export : export data to file")
elif (find(sys.argv[1],"-yplus") >= 0) :
	toto=CSlist()
	for lst in xrange(2,len(sys.argv)) :
		toto.load_listing(sys.argv[lst])
	toto.plot_yplus()
elif (find(sys.argv[1],"-export") >= 0) :
	toto=CSlist()
	toto.load_listing(sys.argv[2])
	toto.export(sys.argv[3])
elif (find(sys.argv[1],"-h") >= 0) :
	print progname+" version "+progversion+" :plot residuals from Code_Saturne listing files"
	sys.exit("Usage: "+progname+" [option] listing1 listing2 ...\n   -yplus : plot yplus min and max\n   -export : export data to file")
else :
	toto=CSlist()
	for lst in xrange(1,len(sys.argv)) :
		toto.load_listing(sys.argv[lst])
	toto.plot_res()
#	toto.plot_ke()