title: Subsetting 2-D latitudes, longitudes and data from a big dataset in Python (IMS snowcover)
body: <p>I have a "simple" objective in mind which turned out to be quite difficult to achieve. I have three 2-Dimensional arrays <code>lat</code>, <code>lon</code> and <code>data</code> all with dimensions (24576, 24576). The first two are latitudes and longitudes coordinates of the variable <code>data</code> which I'm trying to plot over a map. I'm reading all these data from a combination of binary files and text files so that any pre-processing operation is practically impossible and need to be done in the Python script.</p>
  
  <p>Given the dimensionality of the array it is virtually impossible, because of the memory constraints, to directly plot the data even when choosing a basemap projection on a small area over the globe. I've already tried that and getting a memory error when trying to plot with <code>basemap.contourf</code>.</p>
  
  <p>Thus, I need to subset the array BEFORE passing it to the contour function. I've tried many things, but nothing seems to work. My idea was to do something like this</p>
  
  <pre><code>lat_bnds, lon_bnds = [35, 50], [5, 20]
  condition=((lats &gt; lat_bnds[0]) &amp; (lats &lt; lat_bnds[1])) &amp; (lons &gt; lon_bnds[0]) &amp; (lons &lt; lon_bnds[1])
  </code></pre>
  
  <p>Where <code>lats</code> and <code>lons</code> are the 2-D coordinate arrays. This yields a 2-D boolean array which has the same dimension as <code>lats</code> (or equivalently <code>lons</code>) that I can then use to mask the original data but not to subset it. </p>
  
  <p>Using the same condition in <code>numpy.where</code> produces an array shaped <code>(2, 2564856)</code> which I cannot use. I think the problem here is that there are multiple conditions that need to be satisfied in every point of the 2-D array and there is no guarantee that this will lead to a contiguous rectangular submatrix. However for the sake of plotting those data I really need to subset them, or to find another method to elaborate them.</p>
  
  <p>Am I missing something obvious? Is there any other smart way to plot the original data without encountering any error? </p>
  
  <p><strong>Sources of the data:</strong> <a href="http://nsidc.org/data/G02156" rel="nofollow noreferrer">http://nsidc.org/data/G02156</a></p>
  
  <p><strong>Brief script to read the data in assumed that you downloaded the necessary files:</strong></p>
  
  <pre><code>import numpy as np
  import pandas as pd
  from mpl_toolkits.basemap import Basemap # Import the Basemap toolkit
  
  with open('IMS1kmLats.24576x24576x1.double', 'rb') as f:
   data = np.fromfile(f, dtype='&lt;d', count=24576*24576)
   lats = np.reshape(data, [24576, 24576], order='F')
  
  with open('IMS1kmLons.24576x24576x1.double', 'rb') as f:
   data = np.fromfile(f, dtype='&lt;d', count=24576*24576)
   lons = np.reshape(data, [24576, 24576], order='F')
  
  widths=np.full((24576), 1, dtype=int).tolist()
  data=np.array(pd.read_fwf('ims2017312_1km_v1.3.asc', skiprows=30, 
  widths=widths, lineterminator='\n', header=None))
  </code></pre>