!--------------------------------------------------------------------------
   !---------------------------- compute_temperature ----------------------------------
   SELECT CASE(caldyn_thermo)
   CASE(thermo_theta)
      IF(physics_thermo==thermo_fake_moist) THEN
         DO l = ll_begin, ll_end
            !DIR$ SIMD
            DO ij=ij_begin, ij_end
               p_ik = pmid(ij,l)
               theta_ik = temp(ij,l)
               qv = q(ij,l,1) ! water vapor mixing ratio = mv/md
               temp_ik = theta_ik*((p_ik/preff)**kappa)
               temp_ik = temp_ik/(1+0.608*qv)
               temp(ij,l) = temp_ik
            END DO
         END DO
      ELSE
         DO l = ll_begin, ll_end
            !DIR$ SIMD
            DO ij=ij_begin, ij_end
               p_ik = pmid(ij,l)
               theta_ik = temp(ij,l)
               qv = q(ij,l,1) ! water vapor mixing ratio = mv/md
               temp_ik = theta_ik*((p_ik/preff)**kappa)
               temp(ij,l) = temp_ik
            END DO
         END DO
      END IF
   CASE(thermo_entropy)
      IF(physics_thermo==thermo_fake_moist) THEN
         DO l = ll_begin, ll_end
            !DIR$ SIMD
            DO ij=ij_begin, ij_end
               p_ik = pmid(ij,l)
               theta_ik = temp(ij,l)
               qv = q(ij,l,1) ! water vapor mixing ratio = mv/md
               temp_ik = Treff*exp((theta_ik + Rd*log(p_ik/preff))/cpp)
               temp_ik = temp_ik/(1+0.608*qv)
               temp(ij,l) = temp_ik
            END DO
         END DO
      ELSE
         DO l = ll_begin, ll_end
            !DIR$ SIMD
            DO ij=ij_begin, ij_end
               p_ik = pmid(ij,l)
               theta_ik = temp(ij,l)
               qv = q(ij,l,1) ! water vapor mixing ratio = mv/md
               temp_ik = Treff*exp((theta_ik + Rd*log(p_ik/preff))/cpp)
               temp(ij,l) = temp_ik
            END DO
         END DO
      END IF
   CASE(thermo_variable_Cp)
      DO l = ll_begin, ll_end
         !DIR$ SIMD
         DO ij=ij_begin, ij_end
            p_ik = pmid(ij,l)
            theta_ik = temp(ij,l)
            qv = q(ij,l,1) ! water vapor mixing ratio = mv/md
            Cp_ik = nu*( theta_ik + Rd*log(p_ik/preff) )
            temp_ik = Treff* (Cp_ik/cpp)**(1./nu)
            temp(ij,l) = temp_ik
         END DO
      END DO
   CASE(thermo_moist)
      DO l = ll_begin, ll_end
         !DIR$ SIMD
         DO ij=ij_begin, ij_end
            p_ik = pmid(ij,l)
            theta_ik = temp(ij,l)
            qv = q(ij,l,1) ! water vapor mixing ratio = mv/md
            Rmix = Rd+qv*Rv
            chi = ( theta_ik + Rmix*log(p_ik/preff) ) / (cpp + qv*cppv) ! log(T/Treff)
            temp_ik = Treff*exp(chi)
            temp(ij,l) = temp_ik
         END DO
      END DO
   END SELECT
   !---------------------------- compute_temperature ----------------------------------
   !--------------------------------------------------------------------------