diff --git a/src/Emulator.jl b/src/Emulator.jl
index d3fca75ea..8b6f12873 100644
--- a/src/Emulator.jl
+++ b/src/Emulator.jl
@@ -73,6 +73,11 @@ function Emulator(
 
     # For Consistency checks
     input_dim, output_dim = size(input_output_pairs, 1)
+    if isa(obs_noise_cov, UniformScaling)
+        # Cast UniformScaling to Diagonal, since UniformScaling incompatible with 
+        # SVD and standardize()
+        obs_noise_cov = Diagonal(obs_noise_cov, output_dim)
+    end
     if obs_noise_cov !== nothing
         err2 = "obs_noise_cov must be of size ($output_dim, $output_dim), got $(size(obs_noise_cov))"
         size(obs_noise_cov) == (output_dim, output_dim) || throw(ArgumentError(err2))
@@ -240,10 +245,25 @@ function standardize(
     factors::AbstractVector{FT}
 ) where {FT<:AbstractFloat}
     # Case where `output_covs` is a Vector of covariance matrices
+    n = length(factors)
     outputs = outputs ./ factors
     cov_factors = factors .* factors'
     for i in 1:length(output_covs)
-        output_covs[i] = output_covs[i] ./ cov_factors
+        if isa(output_covs[i], Matrix)
+            output_covs[i] = output_covs[i] ./ cov_factors
+        elseif isa(output_covs[i], UniformScaling)
+            # assert all elements of factors are same, otherwise can't cast back to
+            # UniformScaling
+            @assert all(y -> y == first(factors), factors) 
+            output_covs[i] = output_covs[i] / cov_factors[1,1]
+        else
+            # Diagonal, TriDiagonal etc. case
+            # need to do ./ as dense matrix and cast back to original type
+            # https://discourse.julialang.org/t/get-generic-constructor-of-parametric-type/57189/8
+            T = typeof(output_covs[i])
+            cast_T = getfield(parentmodule(T), nameof(T))
+            output_covs[i] = cast_T(Matrix(output_covs[i]) ./ cov_factors)
+        end
     end
     return outputs, output_covs
 end
@@ -268,7 +288,7 @@ dimension). `output_cov` is a Vector of covariance matrices, such as is returned
 function reverse_standardize(
     emulator::Emulator{FT}, 
     outputs::AbstractVecOrMat{FT}, 
-    output_covs::Union{AbstractMatrix{FT}, UniformScaling{FT}, Vector{<:AbstractMatrix{FT}}}
+    output_covs::Union{AbstractMatrix{FT}, Vector{<:AbstractMatrix{FT}}}
 ) where {FT<:AbstractFloat}
     if emulator.standardize_outputs
         return standardize(outputs, output_covs, 1.0 ./ emulator.standardize_outputs_factors)
@@ -297,7 +317,7 @@ in S are sorted in descending order.
 """
 function svd_transform(
     data::AbstractMatrix{FT}, 
-    obs_noise_cov::Union{AbstractMatrix{FT}, UniformScaling{FT}, Nothing}; 
+    obs_noise_cov::Union{AbstractMatrix{FT}, Nothing}; 
     truncate_svd::FT=1.0
 ) where {FT<:AbstractFloat}
     if obs_noise_cov === nothing
@@ -313,7 +333,7 @@ function svd_transform(
         S_cumsum = cumsum(decomp.S) / sum(decomp.S)
         ind = findall(x -> (x > truncate_svd), S_cumsum)
         k = ind[1]
-        n = size(obs_noise_cov)[1]
+        n = size(data)[1]
         println("SVD truncated at k: ", k, "/", n)
         # Apply truncated SVD
 	    transformed_data = sqrt_singular_values_inv[1:k, 1:k] * decomp.Vt[1:k, :] * data
@@ -327,7 +347,7 @@ end
 
 function svd_transform(
     data::AbstractVector{FT}, 
-    obs_noise_cov::Union{AbstractMatrix{FT}, UniformScaling{FT}, Nothing}; 
+    obs_noise_cov::Union{AbstractMatrix{FT}, Nothing}; 
     truncate_svd::FT = 1.0
 ) where {FT<:AbstractFloat}
     # method for 1D data
@@ -337,6 +357,14 @@ function svd_transform(
     return vec(transformed_data), decomposition
 end
 
+function svd_transform(
+    data::AbstractVecOrMat{FT}, obs_noise_cov::UniformScaling{FT}; 
+    truncate_svd::FT = 1.0
+) where {FT<:AbstractFloat}
+    # method for UniformScaling
+    return svd_transform(data, Diagonal(obs_noise_cov, size(data, 1)); truncate_svd=truncate_svd)
+end
+
 """
 svd_reverse_transform_mean_cov(μ, σ2, decomposition::SVD) where {FT}
 
diff --git a/test/Emulator/runtests.jl b/test/Emulator/runtests.jl
index b50d476ee..ef35e5dee 100644
--- a/test/Emulator/runtests.jl
+++ b/test/Emulator/runtests.jl
@@ -8,80 +8,22 @@ using LinearAlgebra
 using CalibrateEmulateSample.Emulators
 using CalibrateEmulateSample.DataContainers
 
-@testset "Emulators" begin
-
-   
-    #build some quick data + noise
-    m=50
-    d=6
-    x = rand(3, m) #R^3
-    y = rand(d, m) #R^5
-    
-    # "noise"
-    μ = zeros(d) 
-    Σ = rand(d,d) 
-    Σ = Σ'*Σ 
-    noise_samples = rand(MvNormal(μ, Σ), m)
-    y += noise_samples
-    
-    iopairs = PairedDataContainer(x,y,data_are_columns=true)
-    @test get_inputs(iopairs) == x
-    @test get_outputs(iopairs) == y
-
-    
-    
-    # [1.] test SVD (2D y version)
-    test_SVD = svd(Σ)
-    transformed_y, decomposition = Emulators.svd_transform(y, Σ, truncate_svd=1.0)
-    @test_throws MethodError Emulators.svd_transform(y, Σ[:,1], truncate_svd=1.0)
-    @test test_SVD.V[:,:] == decomposition.V #(use [:,:] to make it an array)
-    @test test_SVD.Vt == decomposition.Vt
-    @test test_SVD.S == decomposition.S
-    @test size(transformed_y) == size(y)
-
-    # 1D y version
-    transformed_y, decomposition2 = Emulators.svd_transform(y[:, 1], Σ, truncate_svd=1.0)
-    @test_throws MethodError Emulators.svd_transform(y[:,1], Σ[:,1], truncate_svd=1.0)
-    @test size(transformed_y) == size(y[:, 1])
-    
-    # Reverse SVD
-    y_new, y_cov_new = Emulators.svd_reverse_transform_mean_cov(reshape(transformed_y, (d,1)), ones(d,1), decomposition2)
-    @test_throws MethodError Emulators.svd_reverse_transform_mean_cov(transformed_y, ones(d,1), decomposition2)
-    @test_throws MethodError Emulators.svd_reverse_transform_mean_cov(reshape(transformed_y,(d,1)), ones(d), decomposition2)
-    @test size(y_new)[1] == size(y[:, 1])[1]
-    @test y_new ≈ y[:,1]
-    @test y_cov_new[1] ≈ Σ
-    
-    # Truncation
-    transformed_y, trunc_decomposition = Emulators.svd_transform(y[:, 1], Σ, truncate_svd=0.95)
-    trunc_size = size(trunc_decomposition.S)[1]
-    @test test_SVD.S[1:trunc_size] == trunc_decomposition.S
-    @test size(transformed_y)[1] == trunc_size
-    
+#build an unknown type
+struct MLTester <: Emulators.MachineLearningTool end
+
+function constructor_tests(
+    iopairs::PairedDataContainer,
+    Σ::Union{AbstractMatrix{FT}, UniformScaling{FT}, Nothing},
+    norm_factors, decomposition
+) where {FT<:AbstractFloat}
     # [2.] test Normalization
     input_mean = vec(mean(get_inputs(iopairs), dims=2)) #column vector
     sqrt_inv_input_cov = sqrt(inv(Symmetric(cov(get_inputs(iopairs), dims=2))))
-
-    norm_inputs = Emulators.normalize(get_inputs(iopairs),input_mean,sqrt_inv_input_cov)
+    norm_inputs = Emulators.normalize(get_inputs(iopairs), input_mean, sqrt_inv_input_cov)
     @test norm_inputs == sqrt_inv_input_cov * (get_inputs(iopairs) .- input_mean)
 
-    # [3.] test Standardization
-    norm_factors = 10.0
-    norm_factors = fill(norm_factors, size(y[:,1])) # must be size of output dim
-
-    s_y, s_y_cov = Emulators.standardize(get_outputs(iopairs),Σ,norm_factors)
-    @test s_y == get_outputs(iopairs)./norm_factors
-    @test s_y_cov == Σ ./ (norm_factors.*norm_factors')
-    
-    
-    
     # [4.] test emulator preserves the structures
-
-    #build an unknown type
-    struct MLTester <: Emulators.MachineLearningTool end
-    
     mlt = MLTester()
-    
     @test_throws ErrorException emulator = Emulator(
         mlt,
         iopairs,
@@ -92,7 +34,6 @@ using CalibrateEmulateSample.DataContainers
 
     #build a known type, with defaults
     gp = GaussianProcess(GPJL())
-
     emulator = Emulator(
         gp,
         iopairs,
@@ -116,10 +57,9 @@ using CalibrateEmulateSample.DataContainers
         standardize_outputs=false,
         truncate_svd=1.0)
     train_inputs = get_inputs(emulator2.training_pairs)
-    @test norm_inputs == train_inputs
-
-    train_inputs2 = Emulators.normalize(emulator2,get_inputs(iopairs))
-    @test norm_inputs == train_inputs
+    @test train_inputs == norm_inputs
+    train_inputs2 = Emulators.normalize(emulator2, get_inputs(iopairs))
+    @test train_inputs2 == norm_inputs
 
     # reverse standardise
     gp = GaussianProcess(GPJL())
@@ -133,24 +73,106 @@ using CalibrateEmulateSample.DataContainers
         truncate_svd=1.0)
 
     #standardized and decorrelated (sd) data
+    s_y, s_y_cov = Emulators.standardize(get_outputs(iopairs), Σ, norm_factors)
     sd_train_outputs = get_outputs(emulator3.training_pairs)
     sqrt_singular_values_inv = Diagonal(1.0 ./ sqrt.(emulator3.decomposition.S)) 
     decorrelated_s_y = sqrt_singular_values_inv * emulator3.decomposition.Vt * s_y
     @test decorrelated_s_y == sd_train_outputs
+end
 
+@testset "Emulators" begin
+    #build some quick data + noise
+    m=50
+    d=6
+    x = rand(3, m) #R^3
+    y = rand(d, m) #R^5
+    
+    # "noise"
+    μ = zeros(d) 
+    Σ = rand(d,d) 
+    Σ = Σ'*Σ 
+    noise_samples = rand(MvNormal(μ, Σ), m)
+    y += noise_samples
+    
+    iopairs = PairedDataContainer(x,y,data_are_columns=true)
+    @test get_inputs(iopairs) == x
+    @test get_outputs(iopairs) == y
+    
+    # [1.] test SVD (2D y version)
+    test_SVD = svd(Σ)
+    transformed_y, decomposition = Emulators.svd_transform(y, Σ, truncate_svd=1.0)
+    @test_throws MethodError Emulators.svd_transform(y, Σ[:,1], truncate_svd=1.0)
+    @test test_SVD.V[:,:] == decomposition.V #(use [:,:] to make it an array)
+    @test test_SVD.Vt == decomposition.Vt
+    @test test_SVD.S == decomposition.S
+    @test size(transformed_y) == size(y)
 
-    # truncation
-    gp = GaussianProcess(GPJL())
-    emulator4 = Emulator(
-        gp,
-        iopairs,
-        obs_noise_cov=Σ,
-        normalize_inputs=false,
-        standardize_outputs=false,
-        truncate_svd=0.9)
-    trunc_size = size(emulator4.decomposition.S)[1]
-    @test test_SVD.S[1:trunc_size] == emulator4.decomposition.S
-
+    # 1D y version
+    transformed_y, decomposition2 = Emulators.svd_transform(y[:, 1], Σ, truncate_svd=1.0)
+    @test_throws MethodError Emulators.svd_transform(y[:,1], Σ[:,1], truncate_svd=1.0)
+    @test size(transformed_y) == size(y[:, 1])
     
+    # Reverse SVD
+    y_new, y_cov_new = Emulators.svd_reverse_transform_mean_cov(reshape(transformed_y, (d,1)), ones(d,1), decomposition2)
+    @test_throws MethodError Emulators.svd_reverse_transform_mean_cov(transformed_y, ones(d,1), decomposition2)
+    @test_throws MethodError Emulators.svd_reverse_transform_mean_cov(reshape(transformed_y,(d,1)), ones(d), decomposition2)
+    @test size(y_new)[1] == size(y[:, 1])[1]
+    @test y_new ≈ y[:,1]
+    @test y_cov_new[1] ≈ Σ
+    
+    # Truncation
+    transformed_y, trunc_decomposition = Emulators.svd_transform(y[:, 1], Σ, truncate_svd=0.95)
+    trunc_size = size(trunc_decomposition.S)[1]
+    @test test_SVD.S[1:trunc_size] == trunc_decomposition.S
+    @test size(transformed_y)[1] == trunc_size
+    
+    # [3.] test Standardization
+    norm_factors = 10.0
+    norm_factors = fill(norm_factors, size(y[:,1])) # must be size of output dim
 
+    s_y, s_y_cov = Emulators.standardize(get_outputs(iopairs),Σ,norm_factors)
+    @test s_y == get_outputs(iopairs)./norm_factors
+    @test s_y_cov == Σ ./ (norm_factors.*norm_factors')
+    
+    @testset "Emulators - dense Array Σ" begin
+        constructor_tests(iopairs, Σ, norm_factors, decomposition)
+
+        # truncation - only test here, where singular value of Σ differ
+        gp = GaussianProcess(GPJL())
+        emulator4 = Emulator(
+            gp,
+            iopairs,
+            obs_noise_cov=Σ,
+            normalize_inputs=false,
+            standardize_outputs=false,
+            truncate_svd=0.9)
+        trunc_size = size(emulator4.decomposition.S)[1]
+        @test test_SVD.S[1:trunc_size] == emulator4.decomposition.S
+    end
+    
+    @testset "Emulators - Diagonal Σ" begin
+        x = rand(3, m) #R^3
+        y = rand(d, m) #R^5
+        # "noise"
+        Σ = Diagonal(rand(d) .+ 1.0)
+        noise_samples = rand(MvNormal(zeros(d), Σ), m)
+        y += noise_samples
+        iopairs = PairedDataContainer(x, y, data_are_columns=true)
+        _, decomposition = Emulators.svd_transform(y, Σ, truncate_svd=1.0)
+
+        constructor_tests(iopairs, Σ, norm_factors, decomposition)
+    end 
+
+    @testset "Emulators - UniformScaling Σ" begin
+        x = rand(3, m) #R^3
+        y = rand(d, m) #R^5
+        # "noise"
+        Σ = UniformScaling(1.3)
+        noise_samples = rand(MvNormal(zeros(d), Σ), m)
+        y += noise_samples
+        iopairs = PairedDataContainer(x, y, data_are_columns=true)
+        _, decomposition = Emulators.svd_transform(y, Σ, truncate_svd=1.0)
+
+        constructor_tests(iopairs, Σ, norm_factors, decomposition)
+    end 
 end