An extended range of energy stable flux reconstruction schemes, developed using a summation-by-parts approach, is presented on quadrilateral elements for various sets of polynomial bases. For the maximal order bases, a new set of correction functions which result in stable schemes is found. However, for a range of orders it is shown that only a single correction function can be cast as a tensor-product. Subsequently, correction functions are identified using a generalised analytic framework that results in stable schemes for total order and approximate Euclidean order polynomial bases on quadrilaterals -- which have not previously been explored in the context of flux reconstruction. It is shown that the approximate Euclidean order basis can provide similar numerical accuracy as the maximal order basis but with fewer points per element, and thus lower cost.