compute_huber Function

  function compute_huber(this, predicted, expected) result(output)
    !! Compute the huber loss of a model
    implicit none

    ! Arguments
    class(huber_loss_type), intent(in), target :: this
    !! Instance of the huber loss function
    type(array_type), dimension(:,:), intent(inout), target :: predicted
    !! Predicted values
    type(array_type), dimension(size(predicted,1),size(predicted,2)), intent(in) :: &
         expected
    !! Expected values
    type(array_type), pointer :: output
    !! Huber loss

    ! Local variables
    integer :: s, i
    !! Loop indices
    type(array_type), pointer :: ptr
    !! Temporary pointer for calculations

    ptr => predicted(1,1) - expected(1,1)
    output => mean( huber(predicted(1,1) - expected(1,1), this%gamma) )
    if(any(shape(predicted).gt.1))then
       do s = 1, size(predicted,2)
          do i = 1, size(predicted,1)
             if(i.eq.1 .and. s.eq.1) cycle
             if(.not.predicted(i,s)%allocated .or. &
                  .not.expected(i,s)%allocated) cycle
             ptr => predicted(i,s) - expected(i,s)

             output => output + mean( huber(ptr, this%gamma) )
          end do
       end do
    end if

    ! output => merge( &
    !      0.5_real32 * (ptr)**2._real32, &
    !      this%gamma * (abs(ptr) - 0.5_real32 * this%gamma), &
    !      abs(ptr) .le. this%gamma &
    ! )

  end function compute_huber