ftemp_kphio <- function(temp, kphio, kphio_par_a, kphio_par_b){
 out <- kphio * (1.0 + kphio_par_a * (temp - kphio_par_b)^2)
 out <- ifelse(out < 0, 0,
               ifelse(out > 1, 1, 
                      out))
}
library(ggplot2)

ggplot() +
  geom_function(fun = function(x) ftemp_kphio(x, 
                                              kphio = 1, 
                                              kphio_par_a = -0.0015, 
                                              kphio_par_b = 20.0),
                color = "grey") +
  geom_function(fun = function(x) ftemp_kphio(x, 
                                              kphio = 1, 
                                              kphio_par_a = -0.0025, 
                                              kphio_par_b = 20.0)) +
  geom_function(fun = function(x) ftemp_kphio(x, 
                                              kphio = 1, 
                                              kphio_par_a = -0.004, 
                                              kphio_par_b = 20.0),
                color = "grey") +
  geom_hline(yintercept = 0, linetype = "dotted") +
  xlim(1,30) +
  labs(x = "Air temperature (deg C)", y = "Multiplier of the quantum yield")