Dear all, sorry, this might be more a matplotlib question than a kwant question. Below is a minimal code which produces two figures, each showing a different dataset: a current in the scattering region due to a wavefunction in a lead. I would like to combine the two into a single plot, distinguishing the two datasets by e.g. different colormaps (and two different colorbars). Any ideas? Thanks a lot, Tibor ______________________________________ import kwant from math import pi, sqrt from cmath import exp import numpy from matplotlib import pyplot as plt #some parameters t=1 L=30 W=30 EF=0.1 phi=0.01 def hopping(sitei, sitej, phi): xi, yi = sitei.pos xj, yj = sitej.pos return -t*exp(1j*2*pi*phi*2/sqrt(3)*(yj - yi)*(xi + xj)/2) def make_system(): lat = kwant.lattice.general([(sqrt(3)*1/2, 1/2), (0, 1)], [(0, 0), (1/(2*sqrt(3)),1/2)], norbs=1) def central_region(pos): x, y = pos return abs(x) < L/2 and abs(y) < W/2 def lead0_shape(pos): x, y = pos return abs(x) < L/2 #scattering region sys = kwant.Builder() sys[lat.shape(central_region, (0,0))] = -EF sys[lat.neighbors()] = hopping sys.eradicate_dangling() #leads sym = kwant.TranslationalSymmetry(lat.vec((0,1))) lead0 = kwant.Builder(sym) lead0[lat.shape(lead0_shape, (0,0))] = -EF lead0[lat.neighbors()] = hopping lead0.eradicate_dangling() lead1=lead0.reversed() return sys, lead0, lead1 sys, lead0, lead1 = make_system() sys.attach_lead(lead0) sys.attach_lead(lead1) sysf = sys.finalized() wf = kwant.wave_function(sysf, energy=0.0001, params=dict(phi=phi)) J = kwant.operator.Current(sysf) psi1 = wf(0)[0] e_current1 = J(psi1, params=dict(phi=phi)) psi2 = wf(1)[0] e_current2 = J(psi2, params=dict(phi=phi)) kwant.plotter.current(sysf, e_current1, cmap="OrRd", colorbar=True, show=False) ax = plt.gca() plt.colorbar(ax.get_children()[-2], ax=ax) kwant.plotter.current(sysf, e_current2, cmap="Blues", colorbar=True, show=False) ax = plt.gca() plt.colorbar(ax.get_children()[-2], ax=ax) plt.show()