diff --git a/src/display/overloads.jl b/src/display/overloads.jl
index e2c2c98c..21cd192e 100644
--- a/src/display/overloads.jl
+++ b/src/display/overloads.jl
@@ -9,7 +9,8 @@ function show(io::IO, block::InputData)
     Block data:
     """
     dfprintln(io, s)
-    return dfprintln(io, string(block.data) * "\n\n")
+    dfprintln(io, string(block.data) * "\n\n")
+    return 
 end
 
 show(io::IO, data::Vector{InputData}) = map(x -> show(io, x), data)
@@ -25,7 +26,8 @@ function show(io::IO, band::DFBand{T}) where {T<:AbstractFloat}
     string *= """eigvals: $(band.eigvals[1]) -> $(band.eigvals[end])
     extra:   $(band.extra)
     """
-    return dfprintln(io, string)
+    dfprintln(io, string)
+    return 
 end
 
 show(io::IO, bands::Vector{<:DFBand}) = map(x -> show(io, x), bands)
@@ -81,7 +83,6 @@ function show(io::IO, job::DFJob)
         dfprintln(io, crayon"cyan", line, reset)
         dfprintln(io, reset, "(", crayon"green", "scheduled", reset, ", ", crayon"red",
                   "not scheduled", reset, ")")
-        dfprintln(io)
         ln = maximum(length.(string.(name.(is))))
         for (si, i) in enumerate(is)
             n = name(i)
@@ -89,7 +90,8 @@ function show(io::IO, job::DFJob)
             dfprint(io, cr, i == last ? (is_running ? "\t$n <- running\n" : "\t$n <- ran last\n")  : "\t$n\n")
         end
     end
-    return dfprint(io, reset)
+    dfprint(io, reset)
+    return 
 end
 
 function show(io::IO, in::DFCalculation)
@@ -113,14 +115,16 @@ function show(io::IO, in::DFCalculation)
                     "$v\n")
         end
     end
-    return dfprint(io, crayon"reset")
+    dfprint(io, crayon"reset")
+    return 
 end
 
 function show(io::IO, flag_info::QEFlagInfo{T}) where {T}
     df_show_type(io, flag_info)
     dfprintln(io, crayon"yellow", "name : $(flag_info.name)")
     dfprintln(io, crayon"cyan", "description:", crayon"reset")
-    return dfprint(io, "\t" * replace(flag_info.description, "\n" => "\n\t"))
+    dfprint(io, "\t" * replace(flag_info.description, "\n" => "\n\t"))
+    return 
 end
 
 function show(io::IO, flag_info::ElkFlagInfo{T}) where {T}
@@ -133,7 +137,8 @@ function show(io::IO, flag_info::ElkFlagInfo{T}) where {T}
         dfprintln(io, "\t" * "nothing")
     end
     dfprintln(io, crayon"cyan", "description:", crayon"reset")
-    return dfprintln(io, "\t" * replace(flag_info.description, "\n" => "\n\t"))
+    dfprintln(io, "\t" * replace(flag_info.description, "\n" => "\n\t"))
+    return 
 end
 
 function show(io::IO, info::ElkControlBlockInfo)
@@ -143,17 +148,20 @@ function show(io::IO, info::ElkControlBlockInfo)
         dfprintln(io, "\t$(flag.name)")
     end
     dfprintln(io, crayon"cyan", "description:", crayon"reset")
-    return dfprintln(io, "\t" * replace(info.description, "\n" => "\n\t"))
+    dfprintln(io, "\t" * replace(info.description, "\n" => "\n\t"))
+    dfprint(io, crayon"reset")
+    return 
 end
 
 function show(io::IO, el::Element)
     for f in fieldnames(typeof(el))
         dfprintln(io, crayon"red", "$f: ", crayon"reset", "$(getfield(el, f))")
     end
+    dfprint(io, crayon"reset")
 end
 
 function show(io::IO, str::AbstractStructure)
-    dfprintln(io, crayon"cyan", "Structure")
+    dfprintln(io, crayon"cyan", "Structure", crayon"reset")
     dfprintln(io, crayon"red", "    cell parameters:")
     dfprint(io, crayon"reset",
             "\t a = $((str.cell[:,1]...,))\n\t b = $((str.cell[:,2]...,))\n\t c = $((str.cell[:,3]...,))\n")
@@ -162,6 +170,7 @@ function show(io::IO, str::AbstractStructure)
     for a in atoms(str)
         show(io, a)
     end
+    dfprintln(io, crayon"reset")
 end
 
 function show(io::IO, at::AbstractAtom{T,LT}) where {T,LT<:Length{T}}
@@ -183,5 +192,6 @@ function show(io::IO, at::AbstractAtom{T,LT}) where {T,LT<:Length{T}}
             dfprint(io, "$f: $val, ")
         end
     end
-    return dfprintln(io)
+    dfprintln(io, crayon"reset")
+    return 
 end
diff --git a/src/qe/fileio.jl b/src/qe/fileio.jl
index 806fbe8c..c7bdb5da 100644
--- a/src/qe/fileio.jl
+++ b/src/qe/fileio.jl
@@ -42,6 +42,45 @@ function qe_read_output(calculation::DFCalculation{QE}, args...; kwargs...)
     end
 end
 
+function parse_Hubbard_block(f)
+    # Each of these will have n Hubbard typ elements at the end
+    ids = Int[]
+    traces = NamedTuple{(:up, :down, :total), NTuple{3, Float64}}[] 
+    eigvals = (up = Vector{Float64}[], down = Vector{Float64}[])
+    eigvec = (up = Matrix{Float64}[], down=Matrix{Float64}[])
+    occupations = (up = Matrix{Float64}[], down=Matrix{Float64}[])
+    magmoms = Float64[]
+    line = readline(f)
+    cur_spin = :up
+    while strip(line) != "--- exit write_ns ---"
+        line = readline(f)
+        if line[1:4] == "atom"
+            sline = split(line)
+            push!(ids, parse(Int, sline[2]))
+            push!(traces, NamedTuple{(:up, :down, :total)}(parse.(Float64, (sline[end-2], sline[end-1],  sline[end]))))
+            for spin in (:up, :down)
+                readline(f) #should be spin1
+                readline(f)# should be eigvals
+                push!(eigvals[spin], parse.(Float64, split(readline(f))))
+                dim = length(eigvals[spin][1])
+                readline(f) #eigvectors
+                tmat = zeros(dim, dim)
+                for i = 1:dim
+                    tmat[i, :] = parse.(Float64, split(readline(f)))
+                end
+                push!(eigvec[spin], tmat)
+                readline(f) #occupations
+                for i = 1:dim
+                    tmat[i, :] = parse.(Float64, split(readline(f)))
+                end
+                push!(occupations[spin], tmat)
+            end
+            push!(magmoms, parse(Float64, split(readline(f))[end]))
+        end
+    end
+    return [(id = i, trace = t, eigvals = (up = val_up, down = val_down), eigvecs = (up = vec_up, down=vec_down), occupations = (up=occ_up, down=occ_down), magmom = m) for (i, t, val_up,val_down, vec_up,vec_down, occ_up,occ_down, m) in zip(ids, traces, eigvals.up,eigvals.down, eigvec.up,eigvec.down, occupations.up,occupations.down, magmoms)]
+end
+
 """
     qe_read_pw_output(filename::String, T=Float64)
 
@@ -58,7 +97,10 @@ function qe_read_pw_output(filename::String, T = Float64; cleanup = true)
         atsyms       = Symbol[]
         nat          = 0
         while !eof(f)
-            line = readline(f)
+            line = strip(readline(f))
+            if isempty(line)
+                continue
+            end
             #polarization
             if occursin("C/m^2", line)
                 s_line = split(line)
@@ -211,11 +253,29 @@ function qe_read_pw_output(filename::String, T = Float64; cleanup = true)
                 end
                 out[:atomic_positions] = atoms
             elseif occursin("Total force", line)
-                force = parse(T, split(line)[4])
-                if force <= lowest_force
-                    lowest_force      = force
-                    out[:total_force] = force
+                sline = split(line)
+                force = parse(T, sline[4])
+                scf_contrib = parse(T, sline[end])
+                if !haskey(out, :total_force)
+                    out[:total_force] = [force]
+                else
+                    push!(out[:total_force], force)
+                end
+                if !haskey(out, :scf_correction)
+                    out[:scf_correction] = [scf_contrib]
+                else
+                    push!(out[:scf_correction], scf_contrib)
+                end
+            elseif occursin("iteration #", line)
+                sline = split(line)
+                it = length(sline[2]) == 1 ? parse(Int, sline[3]) : parse(Int, sline[2][2:end])
+                if !haskey(out, :scf_iteration)
+                    out[:scf_iteration] = [it]
+                else
+                    push!(out[:scf_iteration], it)
                 end
+                    
+            
             elseif occursin("Magnetic moment per site", line)
                 key = :colin_mag_moments
                 out[key] = T[]
@@ -241,8 +301,10 @@ function qe_read_pw_output(filename::String, T = Float64; cleanup = true)
                     out[key] = [mag]
                 end
             elseif occursin("convergence NOT achieved", line)
-                out[:converged] = false
+                out[:scf_converged] = false
             elseif occursin("convergence has been achieved", line)
+                out[:scf_converged] = true
+            elseif occursin("Begin final coordinates", line)
                 out[:converged] = true
             elseif occursin("atom number", line)
                 if !haskey(out, :magnetization)
@@ -257,16 +319,12 @@ function qe_read_pw_output(filename::String, T = Float64; cleanup = true)
                     out[:magnetization][atom_number] = parse(Vec3{Float64},
                                                              split(readline(f))[3:5])
                 end
-            elseif occursin("Tr[ns(na)]", line)
-                if !haskey(out, :Hubbard_occupation)
-                    out[:Hubbard_occupation] = Float64[]
-                end
-                sline = split(line)
-                id = parse(Int, split(line)[2])
-                if length(out[:Hubbard_occupation]) < id
-                    push!(out[:Hubbard_occupation], parse(Float64, sline[end]))
+            elseif line == "--- enter write_ns ---"
+                
+                if !haskey(out, :Hubbard)
+                    out[:Hubbard] = [parse_Hubbard_block(f)]
                 else
-                    out[:Hubbard_occupation][id] = parse(Float64, sline[end])
+                    push!(out[:Hubbard], parse_Hubbard_block(f))
                 end
                 # Timing info
             elseif occursin("init_run", line)
@@ -384,7 +442,11 @@ function qe_read_pw_output(filename::String, T = Float64; cleanup = true)
                 pop!(out, f, nothing)
             end
         end
-
+        out[:converged] = out[:converged] ? true : haskey(out, :scf_converged) && out[:scf_converged] && !haskey(out, :total_force)
+        if haskey(out, :scf_iteration)
+            out[:n_scf] = length(findall(i -> out[:scf_iteration][i+1] < out[:scf_iteration][i], 1:length(out[:scf_iteration])-1))
+        end
+        pop!(out, :scf_converged, nothing)
         return out
     end
 end
diff --git a/test/testassets/test_job/nscf.in b/test/testassets/test_job/nscf.in
index 98c7826c..21790fa2 100644
--- a/test/testassets/test_job/nscf.in
+++ b/test/testassets/test_job/nscf.in
@@ -3,7 +3,7 @@
   verbosity = 'high'
   calculation = 'nscf'
   outdir = 'localhost/test_job/outputs'
-  pseudo_dir = '/home/ponet/.julia/dev/DFControl/test/testassets/test_job'
+  pseudo_dir = '/home/lponet/.julia/dev/DFControl/test/testassets/test_job'
   restart_mode = 'from_scratch'
 /