_RL theta_Fe_max
theta_Fe_max = theta_Fe_max_lo+
theta_Fe(i,j,k) = theta_Fe_max
& *theta_Fe_max
& /(epsln + alpha_photo2d(i,j,bi,bj)*theta_Fe_max)