cellsOnCell.py

from __future__ import absolute_import, division, print_function # Python 3 Support
import sys

from netCDF4 import Dataset
import numpy as np

if len(sys.argv) != 3: # Argument checking
    print("Usage: python cellsOnCell.py meshFileName cell")
    print("Please supply a mesh file and a cell number of your choosing")
    sys.exit(-1)

meshFileName = sys.argv[1]
cell = int(sys.argv[2])

mesh = Dataset(meshFileName)
nCells = mesh.dimensions['nCells'].size

if cell > nCells:
    print("The requested cell", cell,  "is above nCells", nCells, " of this mesh")
    print("Please supply a cell less then ", nCells)

cellsOnCell = mesh.variables['cellsOnCell'][:]
latCells = mesh.variables['latCell'][:]
lonCells = mesh.variables['lonCell'][:]
nEdgesOnCell = mesh.variables['nEdgesOnCell'][:]
indexToCellID = mesh.variables['indexToCellID'][:]

nEdges = nEdgesOnCell[cell-1] # Use cell - 1 to convert to one base (Python is 0 base) 

print("This requested cell is: ", cell, " and is at ", latCells[cell-1] * (180.0 / np.pi), 
                                                       lonCells[cell-1] * (180.0 / np.pi))
print("And it has: ", nEdges, " neighbors!\n")
print("Its neighbors are...")
print(cellsOnCell[cell-1,:])
for edge in range(nEdges): # A loop over the number of edges/neighbors for this cell 
    neighborCell = cellsOnCell[cell-1, edge]
    neighborLat = latCells[neighborCell-1]
    neighborLon = lonCells[neighborCell-1]

    neighborLat *= (180.0 / np.pi) # Convert to degrees
    neighborLon *= (180.0 / np.pi) # Convert to degrees

    print("\t Neighbor:", edge, " Cell ID: ", neighborCell, "Lat: ", neighborLat, "Lon: ", neighborLon)

output.md

See what cells are neighbors to 1

> python3 cellsOnCell.py x1.2562.grid.nc 1
This requested cell is:  1  and is at  26.56505120270413 -174.9529452257118
And it has:  5  neighbors!

Its neighbors are...
[644 645 646 647 643 643 643 643 643 643]
	 Neighbor: 0  Cell ID:  644 Lat:  27.655334049698872 Lon:  -170.97128347916836
	 Neighbor: 1  Cell ID:  645 Lat:  30.27293040365456 Lon:  -174.95294523008386
	 Neighbor: 2  Cell ID:  646 Lat:  27.655334036761722 Lon:  -178.93460697565868
	 Neighbor: 3  Cell ID:  647 Lat:  23.545808604692223 Lon:  -177.32934203947147
	 Neighbor: 4  Cell ID:  643 Lat:  23.545808610942164 Lon:  -172.57654841257937

Now check one of the edges to see if one of its neighbors is cell id #1:

> python3 cellsOnCell.py x1.2562.grid.nc 644
This requested cell is:  644  and is at  27.655334049698872 -170.97128347916836
And it has:  6  neighbors!

Its neighbors are...
[ 643 1603  163 1606  645    1    1    1    1    1]
	 Neighbor: 0  Cell ID:  643 Lat:  23.545808610942164 Lon:  -172.57654841257937
	 Neighbor: 1  Cell ID:  1603 Lat:  24.426963689883465 Lon:  -168.25096063225894
	 Neighbor: 2  Cell ID:  163 Lat:  28.659673252673446 Lon:  -166.77908290108675
	 Neighbor: 3  Cell ID:  1606 Lat:  31.685750054783977 Lon:  -170.52673583837694
	 Neighbor: 4  Cell ID:  645 Lat:  30.27293040365456 Lon:  -174.95294523008386
	 Neighbor: 5  Cell ID:  1 Lat:  26.56505120270413 Lon:  -174.9529452257118