def plot_ln_abs(gl, *args, **kwargs):
    g = gl.copy()
    for s, b in gl:
        for i in range(len(b.data[0, :, :])):
            for j in range(len(b.data[0, :, :])):
                for n in range(len(b.data[:, 0, 0])):
                    g[s].data[n, i, j] = log(abs(b.data[n, i, j]))
    for s, b in g:
        for i in b.indices:
            for j in b.indices:
                oplot(b[i, j], name = s+'_'+str(i)+str(j), *args, **kwargs)
    plt.gca().set_ylabel('$\\mathrm{ln}\\,\\mathrm{abs}\\,\\tilde{G}(l)$')

def plot_from_archive(archive, function, loops = [-1], indices = [(0, 0)], blocks = ['up'], **kwargs):
    archive = HDFArchive(archive, 'r')
    for l in loops:
        if l < 0:
            ll = archive['results']['n_dmft_loops'] + l
        else:
            ll = l
        for ind in indices:
            for s in blocks:
                f_name = s + '_' + str(ind[0]) + str(ind[1]) + '_it' + str(ll)
                if 'raw' in function: f_name += '_raw'
                f = archive['results'][str(ll)][function]
                if 'iw' in function:
                    if 'RI' in kwargs.keys():
                        if kwargs['RI'] == 'R':
                            oplot(f[s][ind], name = 'Re_' + f_name, **kwargs)
                        if kwargs['RI'] == 'I':
                            oplot(f[s][ind], name = 'Im_' + f_name, **kwargs)
                    else:
                        oplot(f[s][ind], name = 'Re_' + f_name, RI = 'R', **kwargs)
                        oplot(f[s][ind], name = 'Im_' + f_name, RI = 'I', **kwargs)
                elif function == 'g_transf_l':
                    plt.plot(f[s].data[:, ind[0], ind[1]], label = f_name, **kwargs)
                    plt.xlabel('$l_n$')
                    plt.ylabel('$\\tilde{G}(l_n)$')
                else:
                    oplot(f[s][ind], name = f_name, **kwargs)
    y_ax_lab = '$'
    if 'transf' in function: y_ax_lab += '\\tilde{'
    if 'sigma' in function: y_ax_lab += '\\Sigma'
    elif 'g' in function: y_ax_lab += 'G'
    elif 'delta' in function: y_ax_lab += '\\Delta'
    else: y_ax_lab += function
    if 'transf' in function: y_ax_lab += '}'
    y_ax_lab += "("
    if 'iw' in function: y_ax_lab += 'i\\omega_n'
    elif 'tau' in function: y_ax_lab += '\\tau'
    y_ax_lab += ')$'
    plt.gca().set_ylabel(y_ax_lab)
    del archive

from numpy import pi
from matplotlib import pyplot as plt, cm
import sys

n_start = int(sys.argv[1])
n_stop = int(sys.argv[2])
n_step = int(sys.argv[3])
max_w = int(sys.argv[4])
max_a = int(sys.argv[5])

cmaps_pool = [cm.Reds, cm.Blues, cm.Greens, cm.Greys]
for arch in sys.argv[6:len(sys.argv)]:
    print 'loading '+arch+' ...'
    x = CDmft(archive = arch)
    g = x.load('G_sym_iw')
    cmaps = list()
    n_orbs = int(g.n_blocks * .5)
    for i in range(n_orbs):
        cmaps.append(cmaps_pool[i%(len(cmaps_pool))])
    for n in range(n_start, n_stop, n_step):
        for s, b in g:
            g_w = pade(g, s, pade_n_omega_n = n, pade_eta = 10**-2, dos_n_points = 1200, dos_window = (-max_w, max_w), clip_threshold = 0)
            oplot(g_w, RI = 'S', name = s[0]+s[2]+' '+str(n), color = cmaps[int(s[0])](((n - n_start) /float(n_stop - n_start))*.5 + .5))
    filename = 'dos_orb_' + arch[0:-3] + '_' + str(n_start) + str(n_stop) + str(n_step) + '.png'
    plt.gca().set_ylim(0, max_a)
    plt.savefig(filename)
    plt.close()
    print filename + ' ready'
    del x

from pytriqs.gf.local import GfImFreq, SemiCircular

g = GfImFreq(indices = ['eg1','eg2'], beta = 50, n_points = 1000, name = "egBlock")

g['eg1','eg1'] = SemiCircular(half_bandwidth = 1)
g['eg2','eg2'] = SemiCircular(half_bandwidth = 2)

from pytriqs.plot.mpl_interface import oplot, plt
oplot(g, '-o', x_window = (0,10))
plt.ylim(-2,1)

    return 0.7/(z-2.6+0.3*1j) + 0.3/(z+3.4+0.1*1j)

# Semicircle
def GSC(z):
    return 2.0*(z + sqrt(1-z**2)*(log(1-z) - log(-1+z))/pi)

# A superposition of GLorentz(z) and GSC(z) with equal weights
def G(z):
    return 0.5*GLorentz(z) + 0.5*GSC(z)

# Matsubara GF
gm = GfImFreq(indices = [0], beta = beta, name = "gm")
gm << Function(G)
gm.tail.zero()
gm.tail[1] = numpy.array([[1.0]])

# Real frequency BlockGf(reference)
gr = GfReFreq(indices = [0], window = (-5.995, 5.995), n_points = 1200, name = "gr")
gr << Function(G)
gr.tail.zero()
gr.tail[1] = numpy.array([[1.0]])

# Analytic continuation of gm
g_pade = GfReFreq(indices = [0], window = (-5.995, 5.995), n_points = 1200, name = "g_pade")
g_pade.set_from_pade(gm, n_points = L, freq_offset = eta)

# Comparison plot
from pytriqs.plot.mpl_interface import oplot
oplot(gr[0,0], '-o', RI = 'S', name = "Original DOS")
oplot(g_pade[0,0], '-x', RI = 'S', name = "Pade-reconstructed DOS")

#!/usr/bin/env pytriqs
from ClusterDMFT.cdmft import CDmft
from ClusterDMFT.evaluation.analytical_continuation import pade_tr as pade
from pytriqs.gf.local import BlockGf, GfImFreq, GfReFreq
from pytriqs.plot.mpl_interface import oplot
from numpy import pi
from matplotlib import pyplot as plt, cm
import sys

n_start = int(sys.argv[1])
n_stop = int(sys.argv[2])
n_step = int(sys.argv[3])
max_w = int(sys.argv[4])
max_y = int(sys.argv[5])

for arch in sys.argv[6:len(sys.argv)]:
    x = CDmft(archive = arch)
    g = x.load('g_0_c_iw', 0)
    for n in range(n_start, n_stop, n_step):
        g_w = pade(g, pade_n_omega_n = n, pade_eta = 0.001, dos_n_points = 1200, dos_window = (-max_w, max_w), clip_threshold = 0)
        oplot(g_w, RI = 'S', name = str(n), color = cm.jet((n - n_start) /float(n_stop - n_start)))
    filename = 'dosNonInt_' + arch[0:-3] + '_' + str(n_start) + str(n_stop) + str(n_step) + '.png'
    plt.gca().set_ylim(0, max_y)
    plt.savefig(filename)
    plt.close()
    print filename + ' ready'
    del x

def setup_fig():
    axes = plt.gca()
    axes.set_ylabel('$G(\\tau)$')
    axes.legend(loc='lower center',prop={'size':10})

pp = PdfPages('G.pdf')
for plot_objs in objects_to_plot:
    try:
        plt.clf()
        for obj in plot_objs:

            if type(obj) is tuple:
                filename = params.results_file_name(*obj)
                name = 'cthyb'
                if obj[0]: name += '(QN)'
            else:
                filename = 'spinless.' + obj + '.h5'
                name = {'ed':'ED', 'matrix':'Matrix'}[obj]
            
            arch = HDFArchive(filename,'r')
            
            for i1,i2 in product((0,1),(0,1)):
                oplot(arch["tot"][i1,i2], name=name + ",%i%i" % (i1,i2))
                            
        setup_fig()
        pp.savefig(plt.gcf())

    except IOError: pass

pp.close()

from pytriqs.gf.local import *
from pytriqs.plot.mpl_interface import oplot

# A Green's function on the Matsubara axis set to a semicircular
gw = GfImFreq(indices=[1], beta=50)
gw << SemiCircular(half_bandwidth=1)

# Create a Legendre Green's function with 40 coefficients
# initialize it from gw and plot it
gl = GfLegendre(indices=[1], beta=50, n_points=40)
gl << MatsubaraToLegendre(gw)
oplot(gl, "-o")

from pytriqs.gf.local import *
from pytriqs.archive import *
from pytriqs.plot.mpl_interface import oplot

with HDFArchive('slater_five_band.h5','r') as ar:
    # Calculate orbital- and spin-averaged Green's function
    G_tau = ar['G_tau-1']
    g_tau_ave = G_tau['up_0'].copy()
    g_tau_ave.zero()
    for name, g in G_tau: g_tau_ave += g
    g_tau_ave = g_tau_ave/10.
    g_tau_rebin = rebinning_tau(g_tau_ave,1000)
    g_tau_rebin.name = r'$G_{\rm ave}$'
    oplot(g_tau_rebin,linewidth=2,label='')

#!/usr/bin/env pytriqs
from ClusterDMFT.periodization import Periodization
from matplotlib import pyplot as plt
from pytriqs.plot.mpl_interface import oplot
import sys

arch = sys.argv[1]

lat = Periodization(archive = arch)
oplot(_tr_g_lat_pade([lat.get_g_lat_loc()])[0], RI = 'S', name = 'local_DOS')
plt.savefig('plot.png')

from pytriqs.gf.local import *
from pytriqs.archive import HDFArchive
from matplotlib import pyplot as plt
from pytriqs.plot.mpl_interface import oplot

# Read data from archive
ar = HDFArchive('results.h5', 'r')

# Plot imaginary part of the susceptibility on the real axis
oplot(ar['chi_w'][0,0], mode='I', linewidth=0.8, label="$\chi''_0(\\omega)$")
oplot(ar['chi_w'][1,1], mode='I', linewidth=0.8, label="$\chi''_1(\\omega)$")

plt.xlim((-5.0,5.0))
plt.ylim((-1.5,1.5))
plt.ylabel("$\chi(\\omega)$")
plt.legend(loc = "lower right")

from pytriqs.gf.local import *
from pytriqs.archive import *
from pytriqs.plot.mpl_interface import oplot

A = HDFArchive("solution.h5")
oplot(A['Gl']['up'], '-o', x_window=(15,45) )

from pytriqs.gf.local import *
from pytriqs.plot.mpl_interface import oplot

# A Green's function on the Matsubara axis set to a semicircular
gw = GfImFreq(indices=[1], beta=50)
gw << SemiCircular(half_bandwidth=1)

# Create an imaginary-time Green's function and plot it
gt = GfImTime(indices=[1], beta=50)
gt << InverseFourier(gw)
oplot(gt, "-")

    if use_qn: file_name += ".qn"
    file_name += ".h5"

    mkind = lambda spin: (spin,0) if use_blocks else ("tot",spin)

    try:
        arch = HDFArchive(file_name,'r')
        plt.clf()

        name_parts = []
        if use_blocks: name_parts.append('Block')
        if use_qn: name_parts.append('QN')
        name = 'cthyb' + (' (' + ', '.join(name_parts) + ')' if len(name_parts) else '')

        for spin in spin_names:
            bn, i = mkind(spin)
            GF = rebinning_tau(arch['G_tau'][bn],500)
            if use_blocks:
                oplot(GF, name=name + "," + {'up':"$\uparrow\uparrow$",'dn':"$\downarrow\downarrow$"}[spin])
            else:
                i = spin_names.index(i)
                oplot(GF[i,i], name=name + "," + {'up':"$\uparrow\uparrow$",'dn':"$\downarrow\downarrow$"}[spin])
            oplot(ed_arch[spin], name="ED," + {'up':"$\uparrow\uparrow$",'dn':"$\downarrow\downarrow$"}[spin])

        setup_fig()
        pp.savefig(plt.gcf())

    except IOError: pass

pp.close()

from pytriqs.gf.local import *
from pytriqs.archive import HDFArchive
from matplotlib import pyplot as plt
from pytriqs.plot.mpl_interface import oplot

# Read data from archive
ar = HDFArchive('results.h5', 'r')

# Plot input and reconstructed \chi(i\omega_n)
oplot(ar['chi_iw'][0,0],     mode='R', linewidth=0.8, label="$\chi_0(i\\omega_n)$")
oplot(ar['chi_rec_iw'][0,0], mode='R', linewidth=0.8, label="$\chi_\mathrm{0,rec}(i\\omega_n)$")
oplot(ar['chi_iw'][1,1],     mode='R', linewidth=0.8, label="$\chi_1(i\\omega_n)$")
oplot(ar['chi_rec_iw'][1,1], mode='R', linewidth=0.8, label="$\chi_\mathrm{1,rec}(i\\omega_n)$")

plt.xlim((0, 3))
plt.ylabel("$\chi'(i\\omega_n)$")
plt.legend(loc="upper right")

from pytriqs.gf.local import *
from pytriqs.plot.mpl_interface import oplot,plt

# A Green's function on the Matsubara axis set to a semicircular
gw = GfImFreq(indices = [1], beta = 50)
gw <<= SemiCircular(half_bandwidth = 1)

# Create a Legendre Green's function with 40 coefficients
# initialize it from gw and plot it
gl = GfLegendre(indices = [1], beta = 50, n_points = 40)
gl <<= MatsubaraToLegendre(gw)
oplot(gl, '-o')

from pytriqs.gf.local import *
from pytriqs.archive import *
from pytriqs.plot.mpl_interface import oplot

A = HDFArchive("solution.h5")
G_iw = GfImFreq(indices = [0], beta = 50.0)
G_iw << Fourier(A['G_tau']['up'])
oplot(G_iw, '-o', x_window = (0,10))

from pytriqs.plot.mpl_interface import oplot
from pytriqs.gf.local import GfImFreq, Omega, inverse

g = GfImFreq(indices = [0], beta = 300, n_points = 1000, name = "g")
g << inverse( Omega + 0.5 )

# the data we want to fit...
# The green function for omega \in [0,0.2]
X,Y = g.x_data_view (x_window = (0,0.2), flatten_y = True )

from pytriqs.fit import Fit, linear, quadratic

fitl = Fit ( X,Y.imag, linear )
fitq = Fit ( X,Y.imag, quadratic )

oplot (g,     '-o', x_window = (0,5) )     
oplot (fitl , '-x', x_window = (0,0.5) )
oplot (fitq , '-x', x_window = (0,1) )

# a bit more complex, we want to fit with a one fermion level ....
# Cf the definition of linear and quadratic in the lib
one_fermion_level  =  lambda X, a,b   : 1/(a * X *1j  + b),    r"${1}/(%f x + %f)$"    , (1,1)

fit1 = Fit ( X,Y, one_fermion_level )
oplot (fit1 , '-x', x_window = (0,3) )
    

import numpy as np
from pytriqs.gf.local import GfReFreq, SemiCircular

g = GfReFreq(indices = ['eg1', 'eg2'], window = (-5, 5), n_points = 1000, name = "egBlock")

g['eg1','eg1'] = SemiCircular(half_bandwidth = 1)
g['eg2','eg2'] = SemiCircular(half_bandwidth = 2)

from pytriqs.plot.mpl_interface import oplot
oplot(g['eg1','eg1'], '-o', RI = 'S')  # S : spectral function 
oplot(g['eg2','eg2'], '-x', RI = 'S')

#!/usr/bin/env pytriqs
from ClusterDMFT.cdmft import CDmft
from pytriqs.plot.mpl_interface import oplot
import sys
from matplotlib import pyplot as plt

max_freq = int(sys.argv[1])

for n, arch in enumerate(sys.argv[2:]):
    c = CDmft(archive = arch)
    g = c.load('g_c_iw')
    for s, b in g:
        for i in range(len(b.data[0,:,:])):
            for j in range(len(b.data[0,:,:])):
                oplot(b[i, j], x_window = (0, max_freq), RI = 'R', name = 'Re')
                oplot(b[i, j], x_window = (0, max_freq), RI = 'I', name = 'Im')
                if s == 'up': 
                    sout = '\\uparrow'
                elif s == 'down' or s == 'dn': 
                    sout = '\\downarrow'
                else:
                    sout = s
                plt.gca().set_ylabel('$G_{'+sout+'\\,'+str(i)+str(j)+'}^{(c)}(i\\omega_n)$')
                plt.savefig(arch[0:-3]+'_G_'+s+'_'+str(i)+str(j)+'.pdf', dpi = 300)
                plt.close()