classdef Amplifier
    %UNTITLED Summary of this class goes here
    %   Detailed explanation goes here
    
    properties
        type 
        amp_mode
        amplification_db

        nase_mode
        noifig

        saturation_mode
        saturation_power

        fsimu

    end
    
    methods
        function obj = Amplifier(options)
            %UNTITLED Construct an instance of this class
            %   Detailed explanation goes here
            arguments
                options.type
                options.amp_mode
                options.amplification_db

                options.nase_mode = 0;
                options.noifig = 0;

                options.saturation_mode = 0;
                options.saturation_power = 0;

                options.fsimu = []
            end

            obj.type = options.type;
            obj.amp_mode = options.amp_mode;
            obj.amplification_db = options.amplification_db;

            obj.nase_mode = options.nase_mode;
            obj.noifig = options.noifig;

            obj.saturation_mode = options.saturation_mode;
            obj.saturation_power = options.saturation_power;

            obj.fsimu = options.fsimu;


            obj.saturation_power=1/1000*10^(obj.saturation_power/10);


        end
        
        function [y_out,nase] = process(obj,x_in,optional)

            arguments
                obj
                x_in = [];

                optional.nase = 0;
            end

            %calc gain for output power mode, amp mode and saturated mode
            a_lin = obj.calculateGain(x_in);

            y_out = a_lin * x_in;

            if obj.type == "ideal"

                %don't add/ remove noise, but scale it accordingly
                if optional.nase ~= 0
                    nase=state.a^2*optional.nase;
                end

            elseif obj.type == "edfa"

                %calculate ASE-noise
                if optional.nase ~= 0
                    Nase_old = state.gain^2*optional.nase ;

                    h = Constant.Planck;
                    c = Constant.LightSpeed;

                    Nase_new = 0.5* 10^(obj.noifig/10) * h*c/(lambda_T*1e-9) * (a_lin^2-1) ;

                    nase = Nase_old + Nase_new ; 
                end

            end

            pow_in_lin =  mean(abs(x_in.^2)) ; 
            pow_in_dbm = 10*log10(pow_in_lin)+30;
            pow_out_lin = mean(abs(y_out.^2)) ; 
            pow_out_dbm = 10*log10(pow_out_lin)+30;

            if obj.amp_mode == "output_power"
                seemsright = pow_out_dbm == obj.amplification_db;
            elseif obj.amp_mode == "gain"
                seemsright = pow_out_dbm == pow_in_dbm+ obj.amplification_db;
            end

            if ~seemsright 
                warning("Amplifier output not correct, please check the reason");
            end

        end

        function a_lin = calculateGain(obj,xin)

            if obj.amp_mode == "output_power"

                %get linear gain for output power mode
                pow_in = mean(abs(xin.^2)) ; % lin input power
                pow_out = 10^(obj.amplification_db/10 - 3) ; % dBm to lin
                a_lin = sqrt(pow_out/pow_in) ;

            elseif obj.amp_mode == "gain"

                %get linear gain for classic gain mode 
                a_lin=10^(obj.amplification_db/20);
                
            end

            if obj.saturation_mode 
                Pin = sum(mean((abs(xin)).^2, 2)) ;

                gain_power=fzero(inline(['log(G/' num2str(a_lin^2,'%1.16e') ')/(1-G)-log(2)*' num2str(Pin,'%1.16e') '/' ...
                num2str(obj.saturation_power,'%1.16e')],'G'),[1+eps,a_lin^2]);

            end

        end



    end
end