I am interested in truncating a binary MODFLOW CBB file after a certain number of stress periods. Using the Flopy binary file utilities, I have been able to understand how to unpack a binary file and extract values into an array. However, it is unclear if I can use the Flopy utilities directly to create a new binary file. Does Flopy have a utility that would help facilitate writing out a portion of the CBB file to a new binary file? The new CBB file will serve as input to MODPATH, so the formatting of the CBB file would need to remain intact.
With Flopy alone, I feel like I am close to getting the results I need. I'm sure there are Python specific libraries that could be used to get the desired results, but I figured I would see if I can get there with Flopy first.
Here is what I have looked into so far:
import flopy.utils.binaryfile as bf
CBBFile = 'PRE_WT_WP3_PREISS_MidK.cbb'
CBB = bf.CellBudgetFile(CBBFile)
CBB.list_records()
Output from CBB.list_records() for first and last stress periods:
(1, 1, ' STORAGE', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 1, ' CONSTANT HEAD', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 1, 'FLOW RIGHT FACE ', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 1, 'FLOW FRONT FACE ', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 1, 'FLOW LOWER FACE ', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 1, ' DRAINS', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 1, ' RIVER LEAKAGE', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 1, ' HEAD DEP BOUNDS', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 1, ' RECHARGE', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
....
(1, 6, ' STORAGE', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 6, ' CONSTANT HEAD', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 6, 'FLOW RIGHT FACE ', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 6, 'FLOW FRONT FACE ', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 6, 'FLOW LOWER FACE ', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 6, ' DRAINS', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 6, ' RIVER LEAKAGE', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 6, ' HEAD DEP BOUNDS', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
(1, 6, ' RECHARGE', 537, 341, 12, 0, 0.0, 0.0, -1.0, '', '', '', '')
I am looking to keep all of these budget terms in the new binary file, but would like to drop all terms from the last stress period (Stress Period 6 in this example).
I tried to query stress period data into an array and then write this data to a new binary file with the code below. I am running into memory issues with this code. The example CBB file being processed here is considerably smaller (~460 mb) than the CBB I would ultimately like to process (~55 gb). It seems that reading the file to an array may be unnecessary for what I am trying to do.
allRec = CBB.get_data(kstpkper = (0,1))
i = 1
while i < 45:
rec = CBB.get_data(kstpkper = (0,i))
allRec = np.append(allRec, rec)
i += 1
np.save('StrippedCBBFile', allRec)
My other idea was to read the CBB file line-by-line and write out only the data I need. This would require decoding the header line to determine the current stress period to know when to stop writing data. I have not been able to decode the header information to do so.
Thanks!
There are not any specific flopy utilities for rewriting binary budget files. The flopy develop branch has a method (.get_position()) for determining the position of the start of the data (and the header information) for a record in a binary budget file. The .get_position() method will be available in future flopy versions (version > 3.2.10).
>>> import flopy
>>> cobj = flopy.utils.CellBudgetFile('freyberg.gitcbc')
>>> cobj.get_indices(text='CONSTANT HEAD')
array([ 0, 8, 16, ..., 8752, 8760, 8767])
>>> cobj.get_position(8767)
50235424
>>> cobj.get_position(8767, header=True)
50235372
The position data could be used to create a binary file with a slice of the binary budget file using standard python.
fin = open('freyberg.gitcbc', 'rb')
fin.seek(50235372)
length = os.path.getsize(fpth) - 50235372
buffsize = 32
with open('end.cbc', 'wb') as fout:
while length:
chunk = min(buffsize, length)
data = fin.read(chunk)
fout.write(data)
length -= chunk