subroutine control_conv(xin,out)
implicit none
INTEGER, PARAMETER :: iknd=SELECTED_INT_KIND(8)  ! nobody does this
INTEGER, PARAMETER :: DP=SELECTED_REAL_KIND(14)  ! this is double precision
INTEGER, PARAMETER :: SP=SELECTED_REAL_KIND(6)  ! this is single precision  

real(DP) :: Y_hispeed,qcool,CT1,DB,was_hispeed,t_ct_in,t_amb,t_freeze
real(DP) :: Y_closed,Y_fan_closed,T_wb,CT2
real(DP) :: xin(16),out(4),spd1,spd2,Qratio,temp,m,b,t_ct_out
integer :: nspeeds
real(DP),save :: store

    qcool = xin(1)
    CT1 = xin(2) 
    DB = xin(3) 
    T_ct_in = xin(4) 
    T_amb = xin(5) 
    T_freeze = xin(6) 
    t_wb = xin(7) 
    CT2 = xin(8) 
    was_hispeed = xin(9)
    spd1 = xin(10)
    spd2 = xin(11)
    nspeeds = xin(12)
    Qratio = xin(13)
    m = xin(14)
    b = xin(15)
    T_ct_out = xin(16)
    
    Y_CLOSED=0
    Y_FAN_CLOSED=0
    Y_HISPEED=0 

    if(qcool <= 0)then
        Y_closed = 0
        Y_fan_closed = 0
        Y_HISPEED = 0
    elseif(nspeeds == 1)then
        !Set point control
        if(qcool > 0)then
            CT1=CT1-DB/2.*WAS_HISPEED+DB/2.*(1.-WAS_HISPEED)
            IF((T_ct_in >= CT1).AND.(T_AMB > T_FREEZE)) THEN
                Y_CLOSED=1.
                Y_FAN_CLOSED=spd1
                Y_HISPEED=1.
            else
                Y_CLOSED=1.
                Y_FAN_CLOSED=0.1
                Y_HISPEED=0
            endif
        endif           
    elseif(nspeeds == 2)then
        !Differential temperature control
        if(qcool > 0)then
            CT1=CT1-DB/2.*WAS_HISPEED+DB/2.*(1.-WAS_HISPEED)
            IF((T_ct_in >= CT1).AND.(T_AMB > T_FREEZE)) THEN
                Y_CLOSED=1.
                Y_FAN_CLOSED=spd1
                Y_HISPEED=1.
            elseif((T_ct_in > (T_WB+CT2)).AND.(T_AMB > T_FREEZE)) THEN
                Y_CLOSED=1.
                Y_FAN_CLOSED=spd2
                Y_HISPEED=0
            else
                Y_CLOSED=1.
                Y_FAN_CLOSED=0.1
                Y_HISPEED=0
            endif
        endif     
    elseif(nspeeds > 2)then
        !Continuous operation 
        Y_HISPEED=0.
        CT1=CT1-DB/2.*WAS_HISPEED+DB/2.*(1.-WAS_HISPEED)
        if(qcool > 0)then
            Y_closed = 1
            !Y_fan_closed = min(1.,max(0.1,0.736*(T_ct_in - T_ct_out)/(T_ct_in - T_wb) + 0.0399))
            Y_fan_closed = min(1.,max(0.1,0.041*(T_ct_in - T_wb) - 0.05))
            Y_HISPEED = abs(Y_fan_closed == 1)    
        endif       
    endif    
    
    was_hispeed = Y_Hispeed

    out(1) = Y_closed
    out(2) = Y_fan_closed
    out(3) = Y_HISPEED   

end subroutine control_conv