The ability to predict the crystal structure that a molecule will adopt, in advance of the crystallisation experiment or even in advance of synthesis, would have great implications in several areas of materials science. The past two decades have seen great progress in computational methods for crystal structure prediction (CSP), with almost all current methods based on performing a search for the local minima on the high dimensional energy surface representing the energy as a function of the variables that describe a crystal structure. The usual assumption in using these methods is that the global minimum on the potential energy surface corresponds to the most likely observable crystal structure. Focusing on locating the global minimum, many approaches has been developed for CSP, such as Monte Carlo simulated annealing , genetic algorithms and particle swarm optimization. However, higher energy crystal structures are also often observed. This is clear from the prevalence of polymorphism in molecular crystals, where a molecule can adopt more than one crystal structure. Polymorphism is sometimes due to changes in temperature or pressure, which can alter the free energy ordering of structures. These effects can be accounted for in prediction methods by inclusion of entropy and zero-point vibrational contributions to the energy. However, different crystal structures can often be crystallised at the same thermodynamic conditions, sometimes from the same experiment (concomitant polymorphs). It has been estimated, based on a largescale computational study, that nearly 80% pairs of observed polymorphs are monotropic, i.e. their free energies do not cross below their melting temperature. The identification of these metastable polymorphs is important in many applications of CSP, such as polymorph screening of pharmaceuticals, and computer-guided discovery of functional materials, where high energy structures sometimes exhibit the most attractive properties. Thus, for CSP to be predictive of all observable crystal structures of a molecule, the structure search method must not be treated as simply a global energy minimisation problem, but should exhaustively explore the energy landscape for possible structures within the energy range above the global minimum in which observed structures can be located. Therefore, some CSP algorithms, such as lowdiscrepancy, quasi-random sampling place emphasis on exploring the structural landscape as uniformly as possible for all low energy structures. The energy range over which it is important to identify possible crystal structures can be defined by the energy range of observed polymorphism. Most observed polymorphs are separated by only a few kJ mol -1 , although this range can extend above 10 kJ mol -1 in rare cases, or where polymorphs are accessed through desolvation of solvated crystals. Although this defines a narrow energy window for observable structures, the weak interactions between organic molecules mean that large numbers of different structures are often possible within this range, for small organic molecules, a small energy range can include tens or hundreds of putative crystal structures. The importance of the higher energy structures, in addition to the global minimum, creates a tension in designing methods for CSP between efficiency in locating the global energy minimum quickly and time spent exploring the landscape to locate all potentially observable crystal structures. In this work, we describe a hybrid approach, where quasirandom (QR) sampling is used to seed multiple Monte Carlo basin hopping (BH) searches; we refer to the method as QR-BH. The role of quasi-random sampling in QR-BH is to provide a broad sampling of the energy landscape, while basin hopping efficiently locates low energy structures from these starting points. The method is benchmarked on a set of organic molecular crystals and co-crystals to explore its efficiency, how it is influenced by the temperature used in basin hopping and the number of quasi-random seeds vs basin hopping steps used in the search.