#make list of relevant temperatures
temperature=300*Kelvin
#make list of relevant gate voltages
gate_voltage_list=numpy.linspace(-1.0,1.0,9)*Volt
#make list to hold the current calculations
current_list=numpy.zeros(len(gate_voltage_list))
current_list_lin=numpy.zeros(len(gate_voltage_list))
#read the reference transmission spectrum from the netcdf data file
gate_potential_ref = 0.*Volt
bias_potential_ref= 0.2*Volt
transmission_spectrum0 = nlread("gatescan-6-6.nc",
    object_id="trans"+str(gate_potential_ref))[0]
for n in range(len(gate_voltage_list)):
    transmission_spectrum=nlread("gateivscan-6-6.nc",
        object_id="trans"+str(gate_voltage_list[n]))[0]
    current_list[n]=transmission_spectrum.current(
        electrode_temperatures=(temperature,temperature))
    current_list_lin[n]=transmission_spectrum0.current(
        electrode_temperatures=(temperature,temperature),
        electrode_voltages=
            (-0.5*bias_potential_ref+gate_voltage_list[n]-gate_potential_ref,
            0.5*bias_potential_ref+gate_voltage_list[n]-gate_potential_ref)
        )
#plot the current as function of gatevoltage
import pylab
pylab.figure()
pylab.plot(gate_voltage_list,current_list,label="self consistent")
pylab.plot(gate_voltage_list,current_list_lin,label="linear response")
pylab.xlabel("Gate voltage (V)")
pylab.ylabel("Current (A)")
pylab.legend(loc="lower left")
pylab.grid(True)
pylab.show()
