docs: def load_targets(explanatory_rasters):
      """
      Parameters
      ----------
      explanatory_rasters : List of Paths to GDAL rasters containing explanatory variables
  
      Returns
      -------
      expl : Array of explanatory variables
      raster_info : dict of raster info
      """
  
      explanatory_raster_arrays = []
      aff = None
      shape = None
      crs = None
  
      for raster in explanatory_rasters:
          logger.debug(raster)
          with rasterio.open(raster) as src:
              ar = src.read(1)  # TODO band num? 
  
              # Save or check the geotransform
              if not aff:
                  aff = src.affine
              else:
                  assert aff == src.affine
  
              # Save or check the shape
              if not shape:
                  shape = ar.shape
              else:
                  assert shape == ar.shape
  
              # Save or check the geotransform
              if not crs:
                  crs = src.crs
              else:
                  assert crs == src.crs
  
          # Flatten in one dimension
          arf = ar.flatten()
          explanatory_raster_arrays.append(arf)
  
      expl = np.array(explanatory_raster_arrays).T
  
      raster_info = {
          'affine': aff,
          'shape': shape,
          'crs': crs
      }
      return expl, raster_info