Hablas otro idioma?

Friday, 27 May 2011

C++ Example 1: Reading GCM output data (IPCC Data Distribution Centre).

Ejemlo 1 para C++: Lectura de salidas de modelos de circulación global (Datos del IPCC)

Do you know how do climate change studies are carried?

From an engineering perspective, a fundamental element for climate change studies is provided by the outputs from Global Circulation Models, GCM, whose numerical experiments are run for several "climate change scenarios".

Why is it important to learn some programming language? Practitioners (especially young ones) may tend to think that commercial software packages have "the solution". Unfortunately, that is not true, and most of the time we will need to build our own tools. Then, postgraduate students will see the advantages, because those obstacles will enrich their experiences throughout a process where the aim is to understand the basis behind a phenomena.

I am posting a practical example where are described the elements in a C++ program that can be used when dealing with matrices. Why matrices?, because from my perspective, those elements are perhaps the most useful resource in problems where environmental models are involved. You will notice that my programs are not as elaborated as a program written by an IT expert, but they work. Details on the language may be found elsewhere.

At the same time the example shows how to read data obtained from the IPCC Data Distribution Centre.

Glossary:
GCM: Global Circulation Model.
IPCC: Intergovernmental Panel on Climate Change.

EXAMPLE 1: Input matrix: "N" global maps, size "MxP"; output matrix: size "1xN".


Objective: From a global map of near surface temperature, I need to read the data from a single cell (on a certain region of interest).

Issue: The global maps are the output of a Global Circulation Model (GCM), for model runs carried in the period 1960-1999. The outputs of the GCM are provided at monthly resolution (12months * 40years = 480global maps).  I need to retrieve a data from a single cell from each of the 480 global maps.

SOLUTION.

1. Input data format:
2. Using the GRIB converter: Conversion of .grb data into ASCII or bin (thanks to Luminda san´s for the guidance).
  • Using the DOS executable: See readme for full description. In summary, the procedure is:
    • Select output type: e.g., cona (which mean "ASCII output").
    • Type input file: e.g., NCCCSM_20C3M_1_G_tas_1-1560.grb.
    • Type output file: e.g., input.bin.
    • Debug info wanted? yes/no: e.g., no.
    • Wait......the output is generated!!!
  • Format of the output data: See readme for full description. In summary, the data format is:
    • Header  Field Header Field Header Field .....and so on (in rows).
    • Header: 8 columns (the first row of the file).
    • Fields per month: 8 columns; number of rows = (Number of cells on each GCM map matrix) / 8 records of the header = e.g., 128rows*256columns / 8 columns = 4096rows.
    • Total number of rows whole period: e.g., for the period Jan 1870- Dec 1999 (monthly resolution) = (128rows * 256columns) * 12months*130years / 8columns + 12months*130years =  6391320rows.
    • Confirm that the data has been read correctly. I use Fortner´s Transform V3.4, from the USACE Army Geospatial Center, to open the ASCII file.
    • Save the ASCII file into a .bin format, type float (this step is only necessary if you use the code given below).
  • To view, and further edit the arrays you may use the software described in another post.

    3. Reading the data of the region of interest.

    a) Input data.


    // the "//" are used in C to insert comments. They are actually not a part of the code being compiled.
    // the following is the header, which indicates the modules of C that will be used during the compilation.
    #include <stdio.h>
    #include <stdlib.h>


    // the following input variables define the size of the input map. It is shown the case of the CCSM3 GCM.
    // "define" can be used to assign values to permanent variables, which can be called anytime in the programme,
    #define GCMFILA 128  
    #define GCMCOLU 256

    // the following is the format of the input file.
    //1960-1999: (data size of 128*256*12months*40years/8colums)+(header size of 12month*40years)=1966560rows.
    #define FILA 1966560
    #define COLU 8    //data to be read has 8 columns

    //format of the output file
    #define VALUES 480 //outputs at monthly resolution: 12months*40years (1960-1999)
    //mi study area is covered by 1 cell. In case more cells are needed, this value should be increased.
    #define VALUESCOL 4

    // the following are the relative coordinates of the cell I wish to read. For this case I am reading the cell over the Kilimanjaro.
    #define DOWNY 66      
    #define LEFTX 26

    // 128columns*256rows/8 columns+ 1 header row = 4097
    #define FILASMES 4097     

    b) Body of the program.

    main(void)
    {
        // are defined the variables with their corresponding types.
        int i,j,k,m,p;
        float monthlygcm[FILA][COLU];    // this is an auxiliary output file, from which the final results will be read.
        float results[VALUES][VALUESCOL];    // the results are obtained in columns.
        FILE *fp;    // this is called pointer. In this case it is being used to simplify a call to a procedure.

        // this is a common strategy used to open an external file. In this case I am reading a binary file "rb" .
            if((fp=fopen("input.bin","rb"))==NULL){
            printf("cannot open monthlygcm\n");
            exit(1);
        }
       
        // I usually initialize the output matrix.
        for(i=0; i<FILA; i++){
            for(j=0; j<COLU; j++){
                monthlygcm[i][j]=-9999;
            }
        }
       
        //then I call the pointer, which is the input file, and I save it into a matrix which will have all my input data.   
        fread(monthlygcm,sizeof monthlygcm,1,fp);
        fclose(fp);

    //I will read every element of the input matrix "monthlygcm", but I will save only the value of the cell of the region of interest. In this example I am reading the data on the Kilimanjaro, cell (66,26).
    // If you actually plot the GCM output, you will notice that the direction is "upside-down", i.e., the data provided is oriented in the South-North direction, and not in the North-South directions to which we are used to. This means that the actual location of the Kilimanjaro would be the cell (61, 26).


    Please contact me if you need the remaining part of the program. I will only need your name, the use you intent to give to the program, and your compromise to cite the source. The reason to ask that is that I have produced this very simple program during my stay at Tohoku, and I do not want to violate any copyrights. Thanks.

    3 comments:

    1. Hello! My first visit, will visit you again. I am a new follower. Seriously, I thoroughly enjoyed your posts. Congratulations for your work. If you wish to follow back that would be great I'm at http://nelsonsouzza.blogspot.com

      ReplyDelete
    2. Gracias. Poesía. Interesante.

      ReplyDelete
    3. Thank you for stopping by! I really do appreciate your taking the time to leave a comment...I read each and everyone of them. I hope your day is a good one and that you will come back again soon. Take care. Nelson Souzza :)

      ReplyDelete

    Related Posts Plugin for WordPress, Blogger...