code: def qsat(self,temp,pres,parameter):
          
          gamw = (parameter.cl - parameter.cpv) / parameter.rv
          betaw = (parameter.lvtt/parameter.rv) + (gamw * parameter.tt)
          alpw = math.log(parameter.estt) + (betaw /parameter.tt) + (gamw * math.log(parameter.tt))
          work2 = parameter.r/parameter.rv
          foes_lst = [0 for i in range(len(temp))]
          work1_lst = [0 for i in range(len(temp))]
          qsat_lst = [0 for i in range(len(temp))]
  
          for i in range(len(temp)):
            # saturation vapor pressure
            foes_lst[i] = math.exp( alpw - betaw/temp[i] - gamw*math.log(temp[i])  )
            work1_lst[i] = foes_lst[i]/pres[i]
            # saturation humidity
            qsat_lst[i] = work2*work1_lst[i] / (1. + (work2-1.) * work1_lst[i])
  
          return qsat_lst
docs: Calculate (qsat_lst) vector of saturation humidity from:
              temp = vector of element layer temperatures
              pres = pressure (at current timestep).