An Application of CRT Contour Analysis to General Circulation Experiments

Washington W M, O'Lear., B T, Takamine, J, Robertson, D National Center for Atmospheric Research Boulder

September 1968

Bulletin American Meteorological Society

Abstract

The recent development of hemispheric and global general circulation models has produced a need for an efficient method of contour analysis. We will discuss a cathode ray tube (CRT) method of processing contoured data on standard map projections. A detailed explanation of the methods of making black-and-white and color movies is also included.

1. Introduction

The purpose of this paper is to discuss a method of graphical contour analysis of data from numerical experiments on the simulation of the Earth's atmosphere. We have seen within the last few years a development of highly complex hemispheric and global numerical models. Among these are the models of Leith (1965), at the Lawrence Radiation Laboratory, University of California; Mintz (1965), at the University of California, Los Angeles; Smagorinsky, Manabe and Holloway (1965), and Manabe, Smagorinsky and Strickler (1965), at the Geophysical Fluid Dynamics Laboratory, Environmental Science Services Administration; and, more recently, Kasahara and Washington (1967), at the National Center for Atmospheric Research (NCAR). These models have been used successfully to simulate many of the large-scale features of the Earth's atmosphere. We shall discuss a method to process data generated by these models. The details of the method are somewhat unique to NCAR since they are a function of the type of equipment available. However, we expect that the concepts can be extended to almost any type of graphical plotter. For our purposes, we had available a Control Data Corporation (CDC) 6600 computer and a Data Display 80 (dd80) graphical cathode ray tube (CRT) plotter.

In Section 2 we shall discuss some problems associated with analysis of data from general circulation experiments. In Section 3 a method of contouring data onto a cathode ray tube display device is outlined. Section 4 describes various map projections, with the accompanying transformation equations for the cylindrical, polar, and the orthographic projections. In Section 5 we explain the process by which black-and-white and color movies are generated. Finally, in Section 6 we note some future developments needed to improve the processing of such data.

2. Problems in data analysis

Numerical models of the general circulation generate data at grid points in the vertical and horizontal during each of many time steps. This numerical data must be converted to a form which is easily useable to a meteorologist.

The data analysis problem is brought somewhat into perspective when it is shown that for the two-layer NCAR general circulation model (Kasahara and Washington, 1967) there are five data sets of approximately 38,000 numbers generated in 1 min. The plan for future models and computers, as outlined by Carroll and Wetherald (1967) and Kolsky (1966), is to increase the speed of such a calculation by as much as 500 times. It can only be guessed by what factor the data handling problem will be increased when this speed of calculation becomes an actuality.

Fortunately the time step of a normal experiment, 6 min in our case, generates small differences in the variable data set, and it is not necessary to save the data from every time step. The best approach so far has been to save data at some constant interval of model time for the purpose of observing hourly or daily large-scale atmospheric pattern developments. If the data for a given variable or variables are required at every time step, they can be approximated to a suitable degree of accuracy by using an interpolation scheme to fit in needed data between the data sets that were saved. By this method the continuity of computed flow pattern is maintained, and a great deal of time is saved by not writing a data record on magnetic tape at every time step.

The data sets are not in a suitable form for easy analysis since the model generates numbers which are not scaled and adjusted to the types of units most commonly used. In order to see how well the model is accomplishing its task of simulating the atmosphere, it is necessary to observe the output in a form which can be compared to that used for actual atmospheric data. This is done by transforming the data from their primitive form into quantities which can be represented by a graphical contour map. This map is then superimposed on a given projection of the Earth's surface.

3. Methods of contour analysis

Two basic methods are used for two-dimensional contouring. One method is to follow a single contour from beginning to end throughout the entire field and then to do the same with the next contour line until all contours are plotted. This type of analysis is usually used with pen graph-plotters. (See Fawcett, 1962, for application of this method to numerical weather analysis.) A second method is to compute all of the contours within a subregion of the. field and then to move to another subregion until all of the subregions are completed.

Before we discuss both methods, it will be useful to establish some conditions. The first is that the spatial location of a given point is a function of its two-dimensional subscripts. The program is requested to contour a two-dimensional array Z_i,j, where i runs from 1 to m and j from 1 to n, as shown in Fig. 1. The Z_i,j values can be mapped onto x,y coordinates of the CRT scope by x = f(i,j) and y =g(i,j). The f and g are transformation functions. We shall give several examples of transformation functions in Section 4.

The second condition is that linear interpolation is sufficient to find the endpoints of a contour line between two adjacent Z values. The case of Z values being the same as the contour value must be taken into consideration. We always treat equal values the same as values above the line. For a contour line passing between two Z values, one Z value must be above the contour line value (+), the other below (-). Consider the case where Z_i,j is above the line and Z_i,j+1 is below the line in Fig. 2. The wanted endpoint of the contour line of value C is then determined as follows:

where x_e, y_e denote the position of the contour endpoint before transformation to the CRT coordinates.

It is in the determination of which Z values to interpolate that the two contouring methods differ. The simplest logic is in what we call the cell-at-a-time method, which is efficient for a CRT type of contour analysis since a pen does not have to be raised and lowered into above or below the contour line. If all are above or all below, there is no contour line segment at that level in the cell. As shown in Fig. 3, the remaining l4 possibilities fall into one (or two) of the 6 possible line types from which the contour segment can be determined explicitly. The ambiguous cases (Fig. 3, o and p) are arbitrarily assigned, since it is necessary only to preserve consistency of treatment.

While this contouring scheme is extremely simple and therefore fast, it has the drawback of redundant plotter instructions. For example, as shown in Fig. 4, program logic requires that c be calculated to find segment b, even though c was already calculated when segment a was drawn. This is a serious drawback when using a map projection with a complicated transformation formula.

The second contour scheme eliminates the redundant instructions, but at the price of substantially more computer logic. In this scheme, the beginning of each contour line is found and then the line is completely drawn before going on to the next contour. Our scheme is adopted from Dayhoff (1963). Any contour line must be either wholly interior to the array, or start and end on the boundary. Again, points exactly equal to the contour line value are considered to be above the line.

The array is scanned in five parts, as shown in Fig. 5 first the bottom, then the right side, top, left side, and finally the interior. A contour line beginning as shown in Fig. 5 has associated with it a control point and trailing point. The contour is followed by using the control point and searching the surrounding points to see where the contour crosses into a new cell. We have arbitrarily defined the control point as being above the contour line, and we search the surrounding points in a clockwise manner. Fig. 6 shows the surrounding eight points for the interior of the array. Some of the surrounding points are excluded on the boundaries and corners. During the search, each new end point is tested. If it is below the contour value, we proceed until we find a value above. This means that the contour crosses between the last two surrounding points in the search. We then move the control point to the new cell and repeat the procedure. A contour beginning on the boundary must end somewhere on the boundary. Interior contours end when the control point and the search position are identical. This can be verified by first drawing a grid with contours on it and then following the above procedures.

4. Map projections

Three basic projections are usually employed in numerical simulation experiments of the Earth's atmosphere. The cylindrical equal-spaced projection is the most frequently used at NCAR since it displays the entire globe. A polar equidistant projection and an orthographic projection are also available. Sometimes the data over certain regions of the globe, e.g., the tropics, are of particular interest so that special-purpose projections are adapted to provide a more representative piture of the area.

In the cylindrical equal-spaced projection (Fig. 7 ) the parallels and meridians are equally spaced so that the parallels are their correct distance apart. The x₁,y₁ mapping of a point of latitude and longitude X is given by

where a is the mean radius of the Earth. One advantage of the cylindrical equal-spaced projection over other commonly used cylindrical projections, such as the Mercator, is that the entire globe is visible. The scale is preserved at the equator, and distortion increases toward the poles. The poles are greatly exaggerated since the single polar point has been stretched out over the entire longitudinal range. Since the grid data of the NCAR model are at latitude-longitude intersections, and the cylindrical equal-spaced projection is basically in rectangular coordinates, the data from the model are contoured on the projection without transformation.

The polar equidistant projection (Fig. 8) preserves true distance from a pole. As in the cylindrical equal-spaced projection, the parallels in a polar equidistant projection are their correct distance apart. The transformation for the polar projection is given by:

where φ' is colatitude. In this type of projection, distortion increases with distance from the poles. One disadvantage of using the polar projection with a global general circulation model is that two maps are necessary to view the entire globe. However, this type of projection is very useful for polar hemispheric studies. The polar stereographic projection, a modification of the polar projection, is also commonly used.

The orthographic projection (Fig. 9) shows the globe as observed from infinity, and gives the viewer the illusion of perspective and depth. The spacing of parallels and meridians is not uniform (Steers, 1962). The transformation of a point (X, Y, Z) on the surface of a globe to a point (x₃, y₃, z₃) on a plane that passes through the center of the globe will be given. We shall assume that the perpendicular of the plane is centered over a fixed latitude po and longitude M. The transformation equations are:

Verification of the above equations can be made by referring to the schematic of Fig. 10. The plane of projection is defined by x,y coordinates. As shown in the lower part of Fig. 10, we can rotate counterclockwise the x axis of the plane to the X axis of the sphere through the angle λ+90°, and rotate the y axis to the Z axis through the angle 90°-φ. This gives all of the axes coincidence, and yields the above transformation equations (explained by Baldwin, 1963). We should point out that if the transformed point (x_a,y_a,z_a) has a negative z value, then that point is not visible, i.e., it is beyond the horizon. There is a great deal of distortion around the horizon in the orthographic projection. However, this projection can be centered over any given point on the globe, and it is used when a particular area is of special interest.

When regions such as the tropics are of particular interest, special-purpose projections are constructed. In our work, for example, the cylindrical equal-spaced projection has been modified for tropical studies (Fig. 11). The map is divided into two sections to achieve higher resolution on the CRT plotter, which in our case has a square coordinate system.

5. Method of making black-and-white and color movies

The making of computer-generated motion pictures of numerical weather prediction data is relatively new. However, the use of such movies in computer simulation of hydrodynamical flow is not new. Welch (1966), for example, shows the development of water waves in movie form. There are many reasons for making movies of weather data. For example, one can easily digest many interesting features in movies that cannot be seen by looking at maps once or twice a day.

The process by which we make computer-generated movies of weather data has resulted from experimenting with many different combinations of factors. As mentioned in Section 2, we save all fields on a magnetic tape every 4 hr of model time. We then compute time-smoothed values of the quantity which is being displayed using the following formula:

where Z is the smoothed value, Z is the unsmoothed value, t is time, Δt is the time interval (4 hr in our case), α₁ and β₁ are smoothing parameters equal to 0.25 and 2.0, respectively. This type of smoothing heavily clamps motions with a period of 2Δt, motions which would be quite noticeable in a movie since they cause very rapid changes. The time-smoothed values are then interpolated between two time-intervals, say, t and (t+Δt). In our current movies, we interpolate 20 frames. A space-smoothing is finally applied to the interpolated fields to rid them of small-wavelength noise.

where α₂=1/16, β₂=2 and γ₂=4. The above formula has proved quite satisfactory in ridding the field of small-scale features.

Thus far, the process which we have outlined applies both to black-and-white and to color movies. From this point, procedures for generating the two types of film begin to differ. The black-and-white movies are made by causing consecutive frames to be printed over a given time interval of several days. A continuous film strip, which has all the titles and contoured information, is spliced together. The movie is then reduced from 35-mm film generated on the dd80 to 16-mm film for more convenient use.

The color movies are made from black-and-white movies using various techniques. Each part of the movie which is to have a different color is printed on its own black-and-white frame. For example, if we want to make a movie of the sea-level pressure distribution superimposed on a given map projection, the sea-level pressure will be in yellow and the map projection in blue on the final film.

The same methods are used for generating the contoured field information as are used for black-and-white movies. The difference is in how the initial film is made on the dd80. We make alternate frames of black-andwhite film with one frame having just the map projection and continent outline and the alternate frame having the sea-level contours. This process eliminates visible CRT image-drift and provides absolute registration within the frame in the finished film.

A new 35-mm color intermediate negative film is created from this initial film, using a special-effects optical film printer. For our example, the initial film is run through this printer twice, once for each color. The first pass through is with a blue filter, and the map projection and continent outline are exposed in blue. Both films are then rewound, the printer reset to print the alternate frames of sea-level pressure contours, and the blue filter replaced by the yellow. During the second pass, the sea-level pressure contours are exposed on top of the map projection. By processing the color negative film, a master color negative is created. The 35-mm master color negative can then be printed onto 16-mm color film for the finished film. This scheme, which can be used for several colors, involves one difficulty: the process causes the image to turn white when colors overlap. This, however, is not a severe problem if, say, the map projection and the contours do not overlap excessively. A color process that does not have this limitation awaits the further development of color CRT equipment.

Recently we have been able to include clouds in our color movies. Clouds are included by relating their presence and intensity to either the upward vertical velocity or to the moisture content. If the cloud intensity is large within a unit horizontal mesh cell, many plus (+) signs are plotted, but if the cloud intensity is low, few plus signs are plotted. It does not make much difference whether vertical motion or moisture are used as the cloud indicator as long as we allow for cloud-free regions. For example, we assume that clouds do not exist wherever the vertical velocity is less than 0.5 cm sec^-1. By varying this value slightly, we can change the size of the cloudy area. Earlier it was mentioned that the initial film must be run twice through the special-effects optical film plotter in order to obtain a two-color movie. This part of the processing is unchanged. The clouds are superimposed on the intermediate negative film without the use of a color filter; therefore, the color of the clouds is white. In order to avoid creating distinct plus signs on the negative, the lens of the printer is adjusted to make the plus signs out of focus. The result is fuzzy white clouds which show a good correlation to the largescale weather patterns (see cover photograph).

6. Possible future developments

One of the methods of increasing the information content of graphical plotting is to use colors to indicate interactions of different fields. The method described in this article is to superimpose different colors, but, as we mentioned, some difficulty is encountered wherever colors overlap. What we would like to see developed is a CRT device which will plot in different colors without this limitation.

Another valuable innovation would be the use of holography for graphical displays (Pennington, 1968). Holography allows for realistic photography of three-dimensional objects. If its techniques can be developed to reveal the three-dimensional structure of the atmosphere, we will then be able to increase greatly the information content of a single graphical display.

Acknowledgments.

The authors would like to thank the NCAR computer staff for their help in providing the software necessary to develop the CRT programs. We would also like to acknowledge the advice of Dr. Cecil Leith and others at the Lawrence Radiation Laboratory, University of California, where much of the pioneering work on CRT analysis was done. Fred Walden, of NCAR, handled many problems connected with processing film on the dd80 and with making the black-and-white and color films at Color Services, Colorado Springs, Colo. Robert Hite, of the Laboratory of Atmospheric and Biological Sciences at the Goddard Space Flight Center of the National Aeronautics and Space Administration, graciously provided the program which generates the orthographic map projection. Dr. Akira Kasahara, of NCAR, offered valuable advice on the manuscript.

References

Baldwin, H L, Jr, 1963: Orbit display using the SC-4020. Rept. GD/A 63-0459, General Dynamics/Astronautics, San Diego, 38 pp.

Carroll, A B, R T Wetherald, 1967: Application of parallel processing to numerical weather prediction. JACM, 14, 591-614.

Dayhoff, M 0, 1963: A contour-map program for x-ray crystallography. CACM, 6, 620-622.

Fawcett, E B, 1962: Six years of operational numerical weather prediction. J. Appl. Meteor., 1, 318-332.

Kasahara, A, and W M Washington, 1967: NCAR global general circulation model of the atmosphere. Mon. Wea. Rev, 95, 389402.

Kolsky, H G, 1966: Computer requirements in meteorology. IBM Tech. Rept. No. 38.002, 158 pp.

Leith, C E, 1965: Numerical simulation of the Earth's atmosphere. Methods in Computational Physics, Vol. 4, New York, Academic Press, 1-28.

Manabe, S J, Smagorinsky and R F Strickler, 1965: Simulated climatology of a general circulation model with a hydrological cycle. Mon. Wea. Rev., 93, 769-798.

Mintz, Y, 1965: Very long-term global integration of the primitive equations of atmospheric motion. WMO/IUGG Symp. R and D Aspects of Long-Range Forecasting, Boulder, Colo., 1964, WMO Tech. Note No. 66, 141--167.

Pennington, K S, 1968: Advances in holography. Sci. Amer., 218, 40-48.

Smagorinsky, J, Manabe S, Holloway,J L, Jr, 1965: Numerical results from a nine-level circulation model of the atmosphere. Mon. Wea. Rev., 93, 727-768.

Steers, J A, 1962: An Introduction to the Study of Map Projections. London, London University Press, Ltd., 288 pp.

Welch, J E, 1966: simulation of water waves. Datamation, 12, No 11, 41-47.