Loading CRU data and making a shapefile with monthly averages

1035
13
Jump to solution
07-12-2017 04:52 AM
Highlighted
New Contributor II

I have downloaded the 30 year means (1961-1990) of the CRU data from the IPCC website and want to load it into ARCGIS and make shapefiles for each of the datasets (max temp, min temp, precipitation etc). My issue is that the data is formatted to be laid out as east to west in the file from left to right, and north to south going down the file, all in half degree spacing. This is for January and the same layout repeats below for each month in turn.

Is there a way I can load this into ARCGIS and create a shapefile where each row in the attributes table is a single coordinate and each column is as follows: Latitude, Longitude, January, February, March... etc

The data is downloaded from here. When opened with notepad++ the data displays as follows. N.B I can't zoom out enough to show all of the data but if I were to scroll across it would display the world starting from longitude 0.25 right round to 359.75. From the top it starts at latitude of -89.75 and goes down to 89.75 which represents January. Immediately below that, with no break, February follows, then March etc. 

Zoomed in

Zoomed out view of the CRU.dat N.B I can

Reply
0 Kudos
1 Solution

Accepted Solutions
Highlighted
Esri Esteemed Contributor

Hi Suzie Milner , 

Find attached a ZIP with the 12 rasters (ASCII grid format) created with the script posted earlier.

If this answers your question,  could you mark the post with the correct answer as such?

Kind regards, Xander

View solution in original post

13 Replies
Highlighted
Esri Esteemed Contributor

Can you attach a sample of the data to fully understand what you have?

Reply
0 Kudos
Highlighted
New Contributor II

Thanks. I've attached a couple of images and improved the description. Hopefully it's now clearer.

Reply
0 Kudos
Highlighted
Esri Esteemed Contributor

Hi Suzie Milner , thanks for the additional information.

I came across the following thread: CRU 3.2 datasets .dat format (ASCII) not readable in ArcGIS where the suggest to replace the header by a header used by ASCII grid rasters. I tried this with the mean temperature dataset from 1961 - 1990 and changed:

into:

(I also replaced the "-9999"  by " -999")

When I load the data into ArcGIS using WGS1984 as coordinate system an interesting thing occurs:

So the data starts at 0° and the last pixel will be at 360°W... When I change the header to start at 0 instead of -180, the data is matching, but still not what you would want:

The other thing is the range of values:

Since the original headers states "n_months" = 12, does this mean that it needs to be divided by 12?

I'm going to investigate a little more and will post back the results...

Highlighted
MVP Esteemed Contributor

That is the month I presume... If you can get the data into array form, then the process becomes quite simple.

Consider the following view of 'Earth' as 5 rows and 6 columns with each month having a uniform temperature equal to the month 

first

import numpy as np
a = np.zeros((5,6))            # create the world
t = [a+i for i in range(13)]   # give each month a temperature
ts = np.dstack(t)              # produce a years worth of data
m = np.mean(ts, axis=2)        # take the mean for each cell on Earth

m
Out[31]:                       # expected result is as should be
array([[ 6.0,  6.0,  6.0,  6.0,  6.0,  6.0],
       [ 6.0,  6.0,  6.0,  6.0,  6.0,  6.0],
       [ 6.0,  6.0,  6.0,  6.0,  6.0,  6.0],
       [ 6.0,  6.0,  6.0,  6.0,  6.0,  6.0],
       [ 6.0,  6.0,  6.0,  6.0,  6.0,  6.0]])‍‍‍‍‍‍‍‍‍‍‍‍

This can obviously be made transparent to the user... things like bringing the data in as an ascii, then using RasterToNumPyArray to produce the rasters.  But the principle is simple.

It can also be replicated using cell statistics IF they have the spatial analyst license.

Highlighted
Esri Esteemed Contributor

This little script seems to work:

def main():
    in_file = r'C:\GeoNet\IPCC\ctmp6190.dat'
    out_file = r'C:\GeoNet\IPCC\ctmp6190v3.asc'

    with open(in_file, 'r') as f1:
        with open(out_file, 'w') as f2:
            # write new headder
            f2.write('ncols 720\n')
            f2.write('nrows 360\n')
            f2.write('xllcorner -180\n')
            f2.write('yllcorner -90\n')
            f2.write('cellsize 0.5\n')
            f2.write('NODATA_Value -999\n')
            for in_line in f1.readlines():
                if len(in_line) > 100:
                    length = len(in_line)
                    new_line = in_line[length/2:] + in_line[:length/2] + '\n'
                    new_line = new_line.replace('-9999', ' -999')
                    f2.write(new_line)

if __name__ == '__main__':
    main()

The resulting ASC file can be loaded into ArcGIS and displays correctly:

Highlighted
MVP Esteemed Contributor

Xander Bakker‌ deja vu-ish https://community.esri.com/thread/173806 and its preceeding discussion.  Once the data are converted, then it can be processed.  Thought the idea sounded familiar

Highlighted
New Contributor II

Thanks. This enables me to load the data in as a raster file but that loses all of the monthly data. I think it's therefore just showing the January data, being the top matrix of 360 rows. Is there a way I can load it in as something like a shape file that will enable me to keep all 12 month's data, visible in an attribute table?

Reply
0 Kudos
Highlighted
New Contributor II

I noticed that xander_bakker asked if the results need dividing by 12 as it has 12 months data. I would like to retain all 8 months data individually. I think this means I can't use a raster layer but I can't see how to load the data into arcgis otherwise.

Reply
0 Kudos
Highlighted
New Contributor II

Am I right in thinking that to keep the monthly data my only option is to create 12 rater files this way and then create a table from the 12 layers?

Reply
0 Kudos