import math


def Predict_binned_response(spectrum, spectrum_uncertainties, response, response_uncertainties):
    """This function makes a Cherenkov light intensity prediction based on a binned response function



    Parameters

    ----------

    spectrum : array of floats

        Per-bin gamma-ray emission intensities.

    spectrum_uncertainties : array of floats

        Uncertainties in the per-bin gamma-ray emission intensities.
        
    response : array of floats

        Per-bin gamma-ray Cherenkov light response.
        
    response_uncertainties : array of floats

        Uncertanties in the per-bin gamma-ray Chherenkov light response.


    Returns

    -------

    two float values

        the Cherenkov light intensity prediction and the uncertainty in the prediction

    """
    
    prediction = 0.0
    uncertainty = 0.0
    
    for i in range (0,len(spectrum)):
        #Per bin response
        prediction += spectrum[i] * response[i]
        
        uncertainty_sum = 0
        
        #Check that the spectrum or the response has counts, otherwise let the 
        #prediction be 0 and the uncertainty too.
        if spectrum[i] > 0:
            uncertainty_sum += (spectrum[i]*response[i])**2 * (spectrum_uncertainties[i]/spectrum[i])**2
            
        if response[i] > 0:
            uncertainty_sum += (spectrum[i]*response[i])**2 * (response_uncertainties[i]/response[i])**2
        
        #Sum of square of uncertainty
        uncertainty = uncertainty + uncertainty_sum
    
    #And the square root to obtain a RMSE value.
    uncertainty = math.sqrt(uncertainty)
    
    return prediction, uncertainty