program explicit_model_IV_improved
    implicit none
    
    ! --- Definicion de variables ---
    integer :: k, num_puntos
    double precision :: q, k_boltz, beta, alpha
    double precision :: G, G_ref, CR_PTC, T, T_cell, T_ref
    double precision :: Voc_ref, Isc_ref, Jsc_ref, APV
    double precision :: Vt, n_min, n_max, coef, n
    double precision :: Voc, Isc, I0
    double precision :: Vmp, Imp, Pmp, MPPT
    double precision :: Rpv
    
    ! --- Variables del Nuevo Modelo ---
    double precision :: V, I, P_act
    double precision :: paso_I, log_term
    double precision :: V_thermal_n ! Termino combinado (n * Vt)

    ! --- Parámetros Físicos ---
    q = 1.602e-19
    k_boltz = 1.381e-23
    beta  = -0.0023
    alpha =  0.0006

    ! --- Condiciones Ambientales y del Sistema ---
    G = 960.0        ! Irradiancia base
    G_ref = 1361 !1353.0   ! Irradiancia referencia
    CR_PTC = 20!30!16.92   ! Concentracion

    T_cell = 59.4  
    T = T_cell + 273.15
    T_ref = 25.0 + 273.15

    Jsc_ref = 14.8e-3 !10.8e-3 !16.8e-3 
    APV = 12.0
    
    Voc_ref = 2.63 !3.075 
    Isc_ref = Jsc_ref * APV

    ! --- Resistencia Serie (Requerida por el nuevo modelo) ---
    ! Nota: El modelo original no usaba esto explicitamente en el bucle.
    ! Se asume un valor bajo para CPV. Ajustar si se conoce el dato real.
    Rpv = 0.135d0 

    ! --- Cálculo de parámetros dependientes de G y T ---
    Vt = (k_boltz*T)/q
    
    ! Calculo dinamico del factor de idealidad n
    n_max = 2.0
    n_min = 1.0
    coef = 0.0015
    n = n_max + (n_min - n_max)*exp(-coef*G)

    ! Ajuste de Voc e Isc a condiciones de operacion
    Voc = Voc_ref + ((n*k_boltz*T)/q)*log(G/G_ref) + beta*(T - T_ref)
    Isc = Isc_ref * (G/G_ref) * (1.0 + alpha*(T - T_ref)) * CR_PTC

    ! Calculo de corriente de saturacion inversa (I0)
    I0 = Isc / (exp(Voc/(n*Vt)) - 1.0)

    ! Termino agrupado para la ecuacion (n * Vt)
    V_thermal_n = n * Vt

    ! ============================================================
    ! SECCIÓN INTEGRADA: GENERACIÓN DE TABLA Y BÚSQUEDA DE MPP
    ! (Adaptada al método de barrido de Corriente)
    ! ============================================================
    
    ! Configuración del barrido
    ! Aumentamos puntos para mejor resolucion grafica
    num_puntos = 50  
    
    ! IMPORTANTE: Ahora el paso es de CORRIENTE, no de voltaje
    paso_I = Isc / dble(num_puntos)
    
    ! Inicializar variables de búsqueda de MPP
    Pmp = 0.0d0
    Vmp = 0.0d0
    Imp = 0.0d0

    print *, "--------------------------------------------------------"
    print *, " TABLA DE RESULTADOS (Modelo Explicito I -> V)"
    print *, "--------------------------------------------------------"
    print *, " Corriente (I)      Voltaje (V)      Potencia (W)"
    print *, " -------------      -----------      ------------"

    DO k = 0, num_puntos
        ! 1. Definir la Corriente para este paso (Variable Independiente)
        I = dble(k) * paso_I
        
        ! Protección numérica: Evitar I = Isc exacto para que el logaritmo no de error
        ! (log(0) es indefinido). Nos quedamos un epsilon por debajo.
        if (I >= Isc) I = Isc - 1.0e-6

        ! 2. Ecuación del Nuevo Modelo Explícito: V = f(I)
        ! V = n*Vt * ln((Isc - I + I0) / I0) - I*Rs
        log_term = (Isc - I + I0) / I0

        IF (log_term > 0.0d0) THEN
            V = (V_thermal_n * log(log_term)) - (I * Rpv)
        ELSE
            V = 0.0d0
        END IF
        
        ! 3. Corrección física (Voltaje no puede ser negativo)
        IF (V < 0.0d0) THEN
             V = 0.0d0
             P_act = 0.0d0
             ! Imprimimos y salimos del bucle si el voltaje cae a 0 (Cortocircuito)
             WRITE(*, '(F10.4, T20, F10.4, T40, F10.4)') I, V, P_act
             EXIT 
        END IF
        
        ! 4. Calcular Potencia Instantánea
        P_act = V * I
        
        ! 5. Buscar el Máximo (MPP Tracking)
        IF (P_act > Pmp) THEN
            Pmp = P_act
            Vmp = V
            Imp = I
        END IF

        ! 6. Imprimir fila de la tabla
        WRITE(*, '(F10.4, T20, F10.4, T40, F10.4)') I, V, P_act
    END DO
    
    print *, "--------------------------------------------------------"
    print *, ""

    ! ============================================================
    ! FIN DE SECCIÓN INTEGRADA
    ! ============================================================

    ! Calculo final de MPPT del sistema (Escalado x40 segun original)
    MPPT = Pmp * 40.0d0 

    print *, " RESUMEN DE PARAMETROS CALCULADOS:"
    print *, " Voc (cond. op) = ", Voc, " V"
    print *, " Isc (cond. op) = ", Isc, " A"
    print *, " Rpv (asumida)  = ", Rpv, " Ohms"
    print *, " "
    print *, " PUNTO DE MAXIMA POTENCIA (MPP):"
    print *, " Vmp = ", Vmp, " V"
    print *, " Imp = ", Imp, " A"
    print *, " Pmp = ", Pmp, " W"
    print *, " "
    print *, " SISTEMA COMPLETO (40 unidades):"
    print *, " Potencia Total = ", MPPT, " W"
    print *, " "

end program explicit_model_IV_improved