docs: def disvecinflayers(self, x, y, layers, aq=None):
          '''Returns two arrays of size (len(layers),nparam)
          only used in building equations'''
          if aq is None: aq = self.model.aq.find_aquifer_data(x, y)
          qxqy = self.disvecinf(x, y, aq)  # nparam rows, naq cols
          qx = np.sum(qxqy[0,:,np.newaxis,:] * aq.eigvec, 2).T  # Transpose as the first axes needs to be the number of layers
          qy = np.sum(qxqy[1,:,np.newaxis,:] * aq.eigvec, 2).T
          return np.array((qx[layers], qy[layers]))