<filename>src/plotters_simple/PlotWavepackets.py
"""The WaveBlocks Project

Plot homogeneous and inhomogeneous wavepackets and their coefficients
during the time evolution.

@author: <NAME>
@copyright: Copyright (C) 2010, 2011 <NAME>
@license: Modified BSD License
"""

import sys
from numpy import array, conj, abs, angle, squeeze
from matplotlib.pyplot import *

from WaveBlocks import PotentialFactory
from WaveBlocks import HagedornWavepacket
from WaveBlocks import HagedornWavepacketInhomogeneous
from WaveBlocks import IOManager
from WaveBlocks.Plot import plotcf, stemcf

import GraphicsDefaults as GD


def plot_frames_homogeneous(iom, blockid=0, view=None):
    """
    :param iom: An ``IOManager`` instance providing the simulation data.
    """
    parameters = iom.load_parameters()

    # Number of time steps we saved
    timesteps = iom.load_wavepacket_timegrid(blockid=blockid)
    nrtimesteps = timesteps.shape[0]

    # Initialize a Hagedorn wavepacket with the data
    Potential = PotentialFactory().create_potential(parameters)

    # Retrieve simulation data
    grid = iom.load_grid(blockid="global")
    params = iom.load_wavepacket_parameters(blockid=blockid)
    coeffs = iom.load_wavepacket_coefficients(blockid=blockid)

    # A data transformation needed by API specification
    coeffs = [ [ coeffs[i,j,:] for j in xrange(parameters["ncomponents"]) ] for i in xrange(nrtimesteps)]

    # Hack for allowing data blocks with different basis size than the global one
    # todo: remove when we got local parameter sets
    parameters.update_parameters({"basis_size": coeffs[0][0].shape[0]})

    HAWP = HagedornWavepacket(parameters)
    HAWP.set_quadrature(None)

    # Iterate over all timesteps
    for i, step in enumerate(timesteps):
        print(" Plotting frame of timestep "+str(step))

        # Configure the wavepacket and project to the eigenbasis.
        HAWP.set_parameters(params[i])
        HAWP.set_coefficients(coeffs[i])
        HAWP.project_to_eigen(Potential)

        values = HAWP.evaluate_at(grid, prefactor=True)
        coeffi = HAWP.get_coefficients()

        plot_frame(step, parameters, grid, values, coeffi, index=blockid, view=view)

    print(" Plotting frames finished")


def plot_frames_inhomogeneous(iom, blockid=0, view=None):
    """
    :param iom: An ``IOManager`` instance providing the simulation data.
    """
    parameters = iom.load_parameters()

    # Number of time steps we saved
    timesteps = iom.load_inhomogwavepacket_timegrid(blockid=blockid)
    nrtimesteps = timesteps.shape[0]

    # Initialize a Hagedorn wavepacket with the data
    Potential = PotentialFactory().create_potential(parameters)

    # Retrieve simulation data
    grid = iom.load_grid(blockid="global")
    params = iom.load_inhomogwavepacket_parameters(blockid=blockid)
    coeffs = iom.load_inhomogwavepacket_coefficients(blockid=blockid)

    # A data transformation needed by API specification
    params = [ [ params[j][i,:] for j in xrange(parameters["ncomponents"]) ] for i in xrange(nrtimesteps)]
    coeffs = [ [ coeffs[i,j,:] for j in xrange(parameters["ncomponents"]) ] for i in xrange(nrtimesteps)]

    # Hack for allowing data blocks with different basis size than the global one
    # todo: remove when we got local parameter sets
    parameters.update_parameters({"basis_size": coeffs[0][0].shape[0]})

    HAWP = HagedornWavepacketInhomogeneous(parameters)
    HAWP.set_quadrature(None)

    # Iterate over all timesteps
    for i, step in enumerate(timesteps):
        print(" Plotting frame of timestep "+str(step))

        # Configure the wavepacket and project to the eigenbasis.
        HAWP.set_parameters(params[i])
        HAWP.set_coefficients(coeffs[i])
        HAWP.project_to_eigen(Potential)

        values = HAWP.evaluate_at(grid, prefactor=True)
        coeffi = HAWP.get_coefficients()

        plot_frame(step, parameters, grid, values, coeffi, index=blockid, view=view)

    print(" Plotting frames finished")


def plot_frame(step, parameters, grid, values, coeffs, index=0, view=None, imgsize=(12,9)):
    n = parameters["ncomponents"]
    k = array(range(parameters["basis_size"]))

    # Start new plot
    fig = figure(figsize=imgsize)

    for s in xrange(n):
        y = values[s]
        c = squeeze(coeffs[s])

        # Plot the probability densities
        ax1 = fig.add_subplot(n,2,2*s+1)
        ax1.ticklabel_format(style="sci", scilimits=(0,0), axis="y")
        plotcf(grid, angle(y), conj(y)*y)

        if view is not None:
            ax1.set_xlim(view[:2])
            ax1.set_ylim(view[2:])

        ax1.set_xlabel(r"$x$")
        ax1.set_ylabel(r"$\langle\varphi_"+str(s)+r"|\varphi_"+str(s)+r"\rangle$")

        # Plot the coefficients of the Hagedorn wavepacket
        ax2 = fig.add_subplot(n,2,2*s+2)
        ax2.ticklabel_format(style="sci", scilimits=(0,0), axis="y")
        stemcf(k, angle(c), abs(c))

        # axis formatting:
        m = max(abs(c))
        ax2.set_xlim(-1,parameters["basis_size"])
        ax2.set_ylim(-0.1*m, 1.1*m)

        ax2.set_xlabel(r"$k$")
        ax2.set_ylabel(r"$|c|$")

    fig.suptitle(r"Time $"+str(step*parameters["dt"])+r"$")
    fig.savefig("wavepackets_block"+str(index)+"_"+ (7-len(str(step)))*"0"+str(step) +GD.output_format)
    close(fig)


if __name__ == "__main__":
    iom = IOManager()

    # Read file with simulation data
    try:
        iom.open_file(filename=sys.argv[1])
    except IndexError:
        iom.open_file()

    # The axes rectangle that is plotted
    view = [-2.5, 2.5, -0.1, 2.5]

    # Iterate over all blocks and plot their data
    for blockid in iom.get_block_ids():
        print("Plotting frames of data block '"+str(blockid)+"'")
        # See if we have an inhomogeneous wavepacket in the current data block
        if iom.has_inhomogwavepacket(blockid=blockid):
            plot_frames_inhomogeneous(iom, blockid=blockid, view=view)
        # If not, we test for a homogeneous wavepacket next
        elif iom.has_wavepacket(blockid=blockid):
            plot_frames_homogeneous(iom, blockid=blockid, view=view)
        # There is nothing to plot
        else:
            print("Warning: Not plotting any wavepackets in block '"+str(blockid)+"'!")

    iom.finalize()

<filename>src/plotters_simple/PlotWavepackets.py
"""The WaveBlocks Project

Plot homogeneous and inhomogeneous wavepackets and their coefficients
during the time evolution.

@author: <NAME>
@copyright: Copyright (C) 2010, 2011 <NAME>
@license: Modified BSD License
"""

import sys
from numpy import array, conj, abs, angle, squeeze
from matplotlib.pyplot import *

from WaveBlocks import PotentialFactory
from WaveBlocks import HagedornWavepacket
from WaveBlocks import HagedornWavepacketInhomogeneous
from WaveBlocks import IOManager
from WaveBlocks.Plot import plotcf, stemcf

import GraphicsDefaults as GD


def plot_frames_homogeneous(iom, blockid=0, view=None):
    """
    :param iom: An ``IOManager`` instance providing the simulation data.
    """
    parameters = iom.load_parameters()

    # Number of time steps we saved
    timesteps = iom.load_wavepacket_timegrid(blockid=blockid)
    nrtimesteps = timesteps.shape[0]

    # Initialize a Hagedorn wavepacket with the data
    Potential = PotentialFactory().create_potential(parameters)

    # Retrieve simulation data
    grid = iom.load_grid(blockid="global")
    params = iom.load_wavepacket_parameters(blockid=blockid)
    coeffs = iom.load_wavepacket_coefficients(blockid=blockid)

    # A data transformation needed by API specification
    coeffs = [ [ coeffs[i,j,:] for j in xrange(parameters["ncomponents"]) ] for i in xrange(nrtimesteps)]

    # Hack for allowing data blocks with different basis size than the global one
    # todo: remove when we got local parameter sets
    parameters.update_parameters({"basis_size": coeffs[0][0].shape[0]})

    HAWP = HagedornWavepacket(parameters)
    HAWP.set_quadrature(None)

    # Iterate over all timesteps
    for i, step in enumerate(timesteps):
        print(" Plotting frame of timestep "+str(step))

        # Configure the wavepacket and project to the eigenbasis.
        HAWP.set_parameters(params[i])
        HAWP.set_coefficients(coeffs[i])
        HAWP.project_to_eigen(Potential)

        values = HAWP.evaluate_at(grid, prefactor=True)
        coeffi = HAWP.get_coefficients()

        plot_frame(step, parameters, grid, values, coeffi, index=blockid, view=view)

    print(" Plotting frames finished")


def plot_frames_inhomogeneous(iom, blockid=0, view=None):
    """
    :param iom: An ``IOManager`` instance providing the simulation data.
    """
    parameters = iom.load_parameters()

    # Number of time steps we saved
    timesteps = iom.load_inhomogwavepacket_timegrid(blockid=blockid)
    nrtimesteps = timesteps.shape[0]

    # Initialize a Hagedorn wavepacket with the data
    Potential = PotentialFactory().create_potential(parameters)

    # Retrieve simulation data
    grid = iom.load_grid(blockid="global")
    params = iom.load_inhomogwavepacket_parameters(blockid=blockid)
    coeffs = iom.load_inhomogwavepacket_coefficients(blockid=blockid)

    # A data transformation needed by API specification
    params = [ [ params[j][i,:] for j in xrange(parameters["ncomponents"]) ] for i in xrange(nrtimesteps)]
    coeffs = [ [ coeffs[i,j,:] for j in xrange(parameters["ncomponents"]) ] for i in xrange(nrtimesteps)]

    # Hack for allowing data blocks with different basis size than the global one
    # todo: remove when we got local parameter sets
    parameters.update_parameters({"basis_size": coeffs[0][0].shape[0]})

    HAWP = HagedornWavepacketInhomogeneous(parameters)
    HAWP.set_quadrature(None)

    # Iterate over all timesteps
    for i, step in enumerate(timesteps):
        print(" Plotting frame of timestep "+str(step))

        # Configure the wavepacket and project to the eigenbasis.
        HAWP.set_parameters(params[i])
        HAWP.set_coefficients(coeffs[i])
        HAWP.project_to_eigen(Potential)

        values = HAWP.evaluate_at(grid, prefactor=True)
        coeffi = HAWP.get_coefficients()

        plot_frame(step, parameters, grid, values, coeffi, index=blockid, view=view)

    print(" Plotting frames finished")


def plot_frame(step, parameters, grid, values, coeffs, index=0, view=None, imgsize=(12,9)):
    n = parameters["ncomponents"]
    k = array(range(parameters["basis_size"]))

    # Start new plot
    fig = figure(figsize=imgsize)

    for s in xrange(n):
        y = values[s]
        c = squeeze(coeffs[s])

        # Plot the probability densities
        ax1 = fig.add_subplot(n,2,2*s+1)
        ax1.ticklabel_format(style="sci", scilimits=(0,0), axis="y")
        plotcf(grid, angle(y), conj(y)*y)

        if view is not None:
            ax1.set_xlim(view[:2])
            ax1.set_ylim(view[2:])

        ax1.set_xlabel(r"$x$")
        ax1.set_ylabel(r"$\langle\varphi_"+str(s)+r"|\varphi_"+str(s)+r"\rangle$")

        # Plot the coefficients of the Hagedorn wavepacket
        ax2 = fig.add_subplot(n,2,2*s+2)
        ax2.ticklabel_format(style="sci", scilimits=(0,0), axis="y")
        stemcf(k, angle(c), abs(c))

        # axis formatting:
        m = max(abs(c))
        ax2.set_xlim(-1,parameters["basis_size"])
        ax2.set_ylim(-0.1*m, 1.1*m)

        ax2.set_xlabel(r"$k$")
        ax2.set_ylabel(r"$|c|$")

    fig.suptitle(r"Time $"+str(step*parameters["dt"])+r"$")
    fig.savefig("wavepackets_block"+str(index)+"_"+ (7-len(str(step)))*"0"+str(step) +GD.output_format)
    close(fig)


if __name__ == "__main__":
    iom = IOManager()

    # Read file with simulation data
    try:
        iom.open_file(filename=sys.argv[1])
    except IndexError:
        iom.open_file()

    # The axes rectangle that is plotted
    view = [-2.5, 2.5, -0.1, 2.5]

    # Iterate over all blocks and plot their data
    for blockid in iom.get_block_ids():
        print("Plotting frames of data block '"+str(blockid)+"'")
        # See if we have an inhomogeneous wavepacket in the current data block
        if iom.has_inhomogwavepacket(blockid=blockid):
            plot_frames_inhomogeneous(iom, blockid=blockid, view=view)
        # If not, we test for a homogeneous wavepacket next
        elif iom.has_wavepacket(blockid=blockid):
            plot_frames_homogeneous(iom, blockid=blockid, view=view)
        # There is nothing to plot
        else:
            print("Warning: Not plotting any wavepackets in block '"+str(blockid)+"'!")

    iom.finalize()

The WaveBlocks Project  Plot homogeneous and inhomogeneous wavepackets and their coefficients during the time evolution.  @author: <NAME> @copyright: Copyright (C) 2010, 2011 <NAME> @license: Modified BSD License :param iom: An ``IOManager`` instance providing the simulation data. # Number of time steps we saved # Initialize a Hagedorn wavepacket with the data # Retrieve simulation data # A data transformation needed by API specification # Hack for allowing data blocks with different basis size than the global one # todo: remove when we got local parameter sets # Iterate over all timesteps # Configure the wavepacket and project to the eigenbasis. :param iom: An ``IOManager`` instance providing the simulation data. # Number of time steps we saved # Initialize a Hagedorn wavepacket with the data # Retrieve simulation data # A data transformation needed by API specification # Hack for allowing data blocks with different basis size than the global one # todo: remove when we got local parameter sets # Iterate over all timesteps # Configure the wavepacket and project to the eigenbasis. # Start new plot # Plot the probability densities # Plot the coefficients of the Hagedorn wavepacket # axis formatting: # Read file with simulation data # The axes rectangle that is plotted # Iterate over all blocks and plot their data # See if we have an inhomogeneous wavepacket in the current data block # If not, we test for a homogeneous wavepacket next # There is nothing to plot