Comparison of the (a) true and (b) reconstructed topography of the modern Greenland Ice Sheet, using ICESHEET
This post refers to some hints and tricks to run my ice sheet sheet software, which I called ICESHEET. This was published a few weeks ago in Geoscience Model Development, but due to the chaos of moving and starting a new post-doc in Germany, I have not really had the chance to make some further documentation on it. ICESHEET produces perfectly plastic ice sheets under equilibrium conditions. In reality, ice sheets are probably never in equilibrium, but this is probably the simplest approximation of the form of an ice sheet and not break the laws of physics. For a more rigorous scientific description, I suggest reading the Geoscience Model Development paper, which is open access. You can download the software there, or from the downloads page. In fact, even if you got the download from the Geoscience Model Development, I recommend getting the one from the downloads page, because I forgot to include the shapefiles for the Greenland shear stress (whoops).
Some notes on running this. First off, you need to at a minimum have Generic Mapping Tools, as I make use of it to create the grids and do coordinate transformation so that everything is in Cartesian coordinates. The Greenland sample that is included in the download package also requires some manipulation with Netcdf Tools. I fully recommend using the scripts that I included with the Greenland example, though it is possible to use the program without them as long as you have everything in the correct format.
Here are some steps to getting everything to work. Note that I am assuming you are running Linux or some other Unix based system, have Fortran compilers installed already, and are familiar with the command line.
Compile the programs
The base directory contains all the source code for ICESHEET, plus two additional programs that are used in the scripts. If you are running gfortran, you just need to do this:
make icesheet make nearest_int make reduce_dem
If you are running ifort, you have to go into the makefile and uncomment one line that will add a flag that allows the program to read in the binary files. If you decide to put the binaries in another directory (like somewhere in $PATH), you will have to edit the scripts in the later steps to make sure that they can find them.
Basal Shear Stress
The next thing that needs to be done is to set up the basal shear stress files. In the Greenland example, there is a folder called “shear_stress”, everything you need is in there. There are a couple of programs that need to be compiled first, create_ss_grid and convert_grid. Just type
make create_ss_grid make convert_grid
For the shear stress domains file, you must create them yourself and export them as a polygons csv file. This is pretty easy to do if you use QGIS (open the project file). In the shapefile (which should be formatted as a polygon), you should have a field where you enter the shear stress value in Pascals. Further instructions, from the create_shear_stress_file.sh:
# from QGIS: # > Open shapefile # > Layer dropdown menu > save as # > in the options, change "Symbology Export" to "feature symbology", and "Geometry" as "AS_WKT". # > Change the CRS to WGS 84, or else it will not export the points as latitude/longitude, if that is desired
Save the file (call it polygons.csv) in the base shear_stress directory, not in the QGIS directory. If you are making your own domains, make sure there are no gaps between your polygons, or it will revert to the default values (set to 50000 Pa). I’d recommend make sure your shear stress model covers the entire domain you are working on.
When you finished with this, run the create_shear_stress_file.sh script, and it will produce a file that will later be used in ICESHEET.
The key to this is that it produces two files, one containing a grid of the shear stress domains, numbered by the polygon ID in QGIS, and another with the relation between those polygon ID numbers and their corresponding shear stress values. The later file (which is called shear_stress_domain_values.txt) can be edited directly, so you don’t have to edit the polygons in QGIS for each configuration you want. Simply run convert_grid again with the edited shear_stress_domain_values.txt file for each run that you want. For my North American ice sheet complex problem, I basically had a file for each time interval I was interested in.
Creating input grids
After creating the shear stress files as above, the next step is to go into the Greenland directory and run convert_grid.sh. Before doing this, you have to download the IceBridge BedMachine Greenland, Version 2 dataset. Extract it in the base directory with the script. Provided you have the programs in the base directory compiled, the script should create an outline text file (outline.xyz), base topography binary (modern_topo.bin) and basal shear stress binary (ss.bin) that are used when running ICESHEET. It also creates text files that are required by ICESHEET to describe the binary files (elev_parameters.txt and ss_parameters.txt). These have the basic format:
binary_filename min_x max_x min_y max_y grid_spacing
Note that this script must be run in bash shell (if you are running Debian or its offshoots like Ubuntu, the default shell is dash).
Run ICESHEET
The next step is to run “run.sh”. You can change the settings for the program in this file. In particular, you can change the contour interval and spacing between points along the contour. I have it set to 20 m and 20 km, which should take less than a minute to run on any modern computer. The plot2.sh will automatically create the plot files and put them in the plots directory.
Final remarks
For all of the applications I have used ICESHEET on, I simply have made modifications to the scripts given in the source package. The main thing is that if you are starting from something that is in lat/long (the Greenland example is already in Cartesian coordinates). If you want to run something and are having problems, please feel free to email me (evangowan@gmail.com)! I am more than happy to help out.
In Japan, there is this concept called honne and tatemae (本音 and 建前), where you keep your true thoughts to yourself (honne), and project yourself in a more pleasant manner (tatemae). This concept is pretty hard to grasp for a Canadian, where we tend to wallow in our in our (real or imagined) misery. Living in Japan, this is one of the hardest things to get used to. Even simple the simple everyday complaining that I do can be interpreted as a huge grievance.