get_observations_thc Subroutine

    subroutine get_observations_thc(self, obs_values, domain, properties, &
                                    nodal_temperature, nodal_porosity, nodal_pw)
        implicit none
        class(type_output_observation), intent(inout) :: self
        real(real64), intent(out) :: obs_values(:)
        type(type_domain), intent(inout), optional :: domain
        type(type_properties_manager), intent(inout), optional :: properties
        real(real64), intent(in), optional :: nodal_temperature(:)
        real(real64), intent(in), optional :: nodal_porosity(:)
        real(real64), intent(in), optional :: nodal_pw(:)

        type(type_state) :: state
        integer(int32) :: iObs, group_id
        real(real64), allocatable :: original_temperature(:)
        real(real64), allocatable :: original_porosity(:)
        integer(int32) :: istat

        ! Initialize to zero
        obs_values(:) = 0.0d0
        if (.not. present(nodal_temperature)) return
        if (.not. present(nodal_porosity)) return
        if (.not. present(properties)) return
        if (.not. present(nodal_pw)) state%pressure = 101325.0d0

        if (.not. domain%reordering%get_algorithm_name() == "none") then
            allocate (original_temperature, mold=nodal_temperature)
            allocate (original_porosity, mold=nodal_porosity)
            call domain%reordering%to_original_value(nodal_temperature, original_temperature)
            call domain%reordering%to_original_value(nodal_porosity, original_porosity)
        else
            allocate (original_temperature, source=nodal_temperature)
            allocate (original_porosity, source=nodal_porosity)
        end if

        do iObs = 1, self%num_observations
            state%temperature = original_temperature(self%node_ids(iObs))
            state%porosity = nodal_porosity(self%node_ids(iObs))
            group_id = self%elements(iObs)%e%get_group()
            state%water_content = properties%calc_qw(group_id, state)
            if (state%water_content > state%porosity) state%water_content = state%porosity
            if (state%water_content < 0.0d0) state%water_content = 0.0d0
            obs_values(iObs) = properties%calc_thc(group_id, state)
        end do

        deallocate (original_temperature)
        deallocate (original_porosity)

    end subroutine get_observations_thc