Using Areas
Previous Chapter Tutorial Home Next Chapter
- Ar 1. What is Areas?
- Ar 1.1 Table of Areas user entry points
- Initialization and support routines
- Action routines
- Query routine
- Parameter setting and retrieving routines
- Debugging routine
- Ar 1.2 Table of Areas parameters
- Ar 1.3 Areas definitions: Edges
- Ar 1.4 Areas definitions: Groups of edges
- Ar 1.5 Areas definitions: Left and right
- Ar 1.6 Areas definitions: Areas, area identifiers, and area maps
- Ar 1.7 Areas definitions: Masks and areas defined by more than one edge group
- Ar 1.8 What calls do I need to get my Areas plot?
- Ar 1.9 Areas parameters: What they do and how to use them
- Ar 2. Initializing Areas
- Ar 2.1 Set up area map
- Ar 2.2 Adding edges to an area map
- Ar 2.3 Preprocess area map: Overview
- Ar 2.3.1 Preprocess area map: Breaking up edges
- Ar 2.3.2 Preprocess area map: Removing coincident and dangling edges
- Ar 2.3.3 Preprocess area map: Reconciling area information
- Ar 3. Producing results using Areas
- Ar 3.1 Obtaining area ids and group ids
- Ar 3.2 Masking and drawing lines through an area map
- Ar 3.3 Writing a masking or line-drawing routine
- Ar 3.4 Filling areas
- Ar 3.5 Writing an area-fill routine
- Ar 4. Debugging Areas
- Ar 4.1 Debugging options: Setting parameters
- Ar 4.2 Debugging options: Drawing parts of the area map
- Ar 5. Areas parameter descriptions
Areas is designed as a support utility for other utilities. For this reason, most users need Areas only when using Conpack, Ezmap, or Vectors. However, Areas has a broad range of applications beyond the uses discussed in this chapter.
Used with a utility like Ezmap or Conpack that creates lines dividing the plotter frame into regions, Areas allows you to retrieve (or recover) outline and other information about each region. Also, Areas allows you to overlay two or more groups or sets of lines, and to recover the regions defined by the union of the groups as well as telling you where each region is located relative to each group of lines. These features allow you to do things like mask contour lines or latitude/longitude lines over the ocean, and color-fill or pattern-fill regions.
The Areas utility has nine user entry points. This table organizes the routines according to their functions.
- ARINAM
- Initializes an area map. ARINAM must be called before using any
other Areas routine that references that area map. See modules Ar 2.1, Mp
4.3, and Cp
5.3.
- AREDAM
- Adds edges or area boundaries to an area map. AREDAM must be
called for each edge to be included in an area map. See module Ar 2.2.
- ARPRAM
- Processes an area map, removing nonclosed edges and cleaning up
area identifiers, among other things. ARPRAM is called by the action
and query subroutines if you don't call it explicitly. See modules Ar 2.3 through Ar 2.3.3.
- ARDRLN
- Breaks a given polyline into pieces, each of which lies in a
single area in an area map. A specified user routine is then used to
decide what to do with each piece. See module Ar
3.2.
- ARSCAM
- Scans an area map, extracting the definitions of all areas formed
by the map, and allows you to process each area separately. See
modules Ar 3.4, Mp 4.7, and Cp 5.7.
- ARGTAI
- Given a point in an area in an area map, ARGTAI returns group and area information about that area. See module Ar 3.1.
- ARGETI
- Retrieves an Areas integer parameter value. See module Ar 1.9.
- ARSETI
- Sets an Areas integer parameter value. See module Ar 1.9.
- ARDBPA
- Produces a picture of that part of the area map that belongs to a specified group. See module Ar 4.2.
This table is a quick reference guide to all Areas parameters. More
complete descriptions appear in section "Ar 5.
Areas parameter descriptions."
The behavior of a typical routine in an NCAR Graphics utility is
sometimes determined entirely by the routine's arguments, but
frequently it is also affected by the value of one or more of the
utility's "parameters." A parameter is a variable that
controls the behavior of a utility; parameters are accessed via
parameter-access routines that can set or retrieve the parameter
value.
Instructions for setting and retrieving Areas parameters are
provided in module "Ar 1.9 Areas parameters: What
they do and how to use them."
-------------------------------------------------------------------
Parameter Brief description Fortran type Examples Module
-------------------------------------------------------------------
AT Arithmetic Type flag Integer
DB DeBug plots flag Integer cardb1 Ar 4.1
DC Debug Colors index Integer
DI Direction Indicator flag Integer
LC Largest Coordinate Integer
-------------------------------------------------------------------
The Areas definitions modules Ar 1.3 through
Ar 1.7 provide the background you need to
understand the discussions in subsequent Areas modules. Because
technically complete definitions involve complex mathematics that are
not needed by most users, the definitions provided here are intended
to give an intuitive understanding of the material.
Three related concepts are very basic to Areas: edges, groups of edges, and group identifiers.
Figure 1 Figure 2
Figure 3 Figure 4
An edge is a set of line segments connected so that a pen tracing the edge never needs to be lifted. Said another way, points specifying an edge must be sequentially ordered, as in Figure 1. Edges can be made of more than one line segment.
When one or more edges cross, Areas adds a vertex at each intersection so that regions such as A and B in Figure 2 can be addressed separately. When Figure 2 is drawn, it is defined by line segments (1,2), (2,3), (3,4), (4,1), but when the edges are processed by Areas, the point Z is added, and the edge segments defining the figure become (1,Z), (Z,2), (2,3), (3,Z), (Z,4), (4,1). These edge segments are among the building blocks of the Areas utility.
It is possible to define the same set of edge segments in many different ways. Looking at Figure 3, you can define the square as a single edge (1, 2, 3, 4, 1), or as two or more edges, like (1, 2, 3) and (3, 4, 1). In both cases, the square is defined by a set of edges: an edge group.
The most useful edge groups have no vertexes from which only one line segment emanates. In other words, useful edge groups won't have dangling line segments. This property implies that the plane is divided into discrete areas with well-defined boundaries. For example, the pentagon in Figure 4 divides the plotter frame into the area inside the pentagon and the area outside the pentagon. Notice that the line segment between points 1 and 2 in Figure 4 does have a vertex with precisely one line segment emanating from it. By default, Areas will remove dangling edges like this from the area map.
An area map is a linked list of edge segments, from which area and
group information can be extracted. Area maps are discussed in module
"Ar 1.6 Areas definitions: Areas, area
identifiers, and area maps."
Ezmap produces dangling line segments. For instance, some small islands are crescent-shaped edges, so Areas will remove them from the area map.
A group of edges can be just about any set of edges that you
want to group together. The most useful groups of edges are boundary
lines from the Ezmap utility, contours from Conpack, and lines
defining vertical strips. By defining more than one group of edges in
an area map, it is possible to access areas defined by each group as
well as areas defined by the union of the groups.
Note: Areas definition modules Ar 1.3
through Ar 1.7 provide the background needed
for the discussions in subsequent Areas modules. Because technically
complete definitions involve complex mathematics that are not needed
by most users, the definitions provided here are intended to give an
intuitive understanding of the material.
The graphic in this module shows three sets or groups of edges. Areas
allows you to individually assign any edge to an edge group. Groups
are differentiated using a positive integer value called a group
identifier. The module "Ar 2.2 Adding
edges to an area map" discusses the details of assigning
edges and group identifiers.
Although you can assign any identifier to a specific group with Areas, it is helpful to know the group identifiers for the groups that are created with the Ezmap and/or Conpack utilities. The following group identifiers are the defaults used with Ezmap and Conpack:
- 1 Regions defined by the Ezmap outline set
- 2 Vertical strips defined by Ezmap
- 3 Contour levels defined by Conpack
- 4 Vertical strips defined by Conpack
Vertical strips are strips dividing an Ezmap or Conpack plot into columns so that it will plot correctly on limited hardware devices, such as older terminals.
Using the graphic in this module, you could overlay the contour map over the geographic map; the map outlines would all be in group 1, and the contour lines would all be in group 3. By default, Ezmap produces a single vertical strip for the entire plot, so you would also have a group 2 defined that contains the perimeter of the geographic map. If you needed to overlay Conpack vertical strips as well, the Conpack vertical strips are placed in group 4 by default.
Each edge segment has a left side and a right side.
Note: Definitions in modules Ar 1.3
through Ar 1.7 provide the background needed
for the discussions in subsequent Areas modules. Because technically
complete definitions involve complex mathematics that are not needed
by most users, the definitions provided here are intended to give an
intuitive understanding of the material.
Figure 1a Figure 1b
Figure 2 Figure 3
Every edge segment has a left side and a right side that depends on
the order of the points in that edge. The best way to understand which
side is left and which is right is to imagine yourself standing on the
first point on an edge, facing the second. Everything to your left is
left of the edge segment, and everything to your right is right of the
edge segment. Notice how in Figure 1, the area to the right of the
edge depends on the order of the points in the edge.
Figure 2 shows a more complex edge. In this case, the area inside the
star is to the left of the edge, and the area outside the star is to
the right of the edge. Figure 3 illustrates an edge that crosses
itself. Areas breaks the edge down into edge segments, and then in
sequential order, determines what is right of each edge segment and
what is left of each edge segment. The preprocess area map modules Ar 2.3 through Ar 2.3.3 discuss
how area information is reconciled for edges like the one in Figure
3.
Conpack always traces contour lines in such a way that the function
being contoured slopes downhill to the left and uphill to the right.
For this reason, areas where the function value is below the contour
line are to the left and areas where the function value is above the
contour line are to the right. Since this is always the case, we refer
to areas above or below a contour line rather than to the left or
right of a contour line.
This module clarifies how Areas defines the term area; this
module also defines the term area map.
Note: Definitions in modules Ar 1.3--Ar 1.7 offer an intuitive approach to
understanding the material.
Figure 1 Figure 2a
Figure 2b Figure 3a
Figure 3b Figure 4
An area defined by a set E of edge segments is a set S of
points in the plane such that 1) if A and B are both members of S,
there exists a continuous curve that connects A and B and that nowhere
crosses any of the edge segments in E, and 2) if A is a member of S
and there does exist a continuous curve that connects A to some
other point B without crossing any of the edge segments in E, then B
is also a member of S.
Let's look at some cases illustrating what this definition means.
First, Areas defines a perimeter for the plotter frame so areas
that would be otherwise unbounded can be identified and treated
correctly by limited software and hardware. This perimeter is drawn
around the limits of your plotter frame, but you can make any edge
group inside your plotter frame into a perimeter. Figure 1 shows the
default perimeter created by Areas.
First let's assume that we have two points, A and B, with an edge
segment joining them. As shown in Figure 2a, the area defined by A and
B is the entire plotter frame. When both A and B lie on different
edges of the perimeter (Figure 2b), the edge segment (A, B) defines
two areas, one to the left of the edge and one to the right of the
edge.
A more complex case is that of a closed edge, or a loop (Figure 3a).
In this case, what is outside the loop is an area, and what is inside
the loop is an area. By adding one or more loops inside or outside the
original, we now have areas within areas, or alternatively, holes
within areas (Figure 3b).
Another common edge group that defines areas is a set of intersecting
edges, for example, a map of the US state boundaries. As shown in
Figure 4, you can think of an area as a region that is defined by
starting at an endpoint X on an edge segment, going to the end of that
edge segment, choosing the leftmost connected edge segment, going to
the end of it, choosing the leftmost connected segment, and so on,
until you return the first edge segment (Figure 4). To define another
area, start at an edge segment and choose only rightmost edges.
Each edge segment has a left area identifier and a right area
identifier associated with it. As an area is defined, Areas keeps
track of the area identifiers associated with each edge segment; if
these don't all match, Areas decides which one has precedence. This
process is discussed in detail in the preprocess area map modules 2.3 through 2.3.3.
Generally, negative area identifiers are all treated as if they were -1; this denotes areas outside the perimeter. Zero area identifiers are used to denote areas that have not yet had an area identifier assigned, or areas where you don't know what area identifier has been assigned by some other edge segment in the area map. Positive area identifiers should be used to denote any area that you may want to address. The numbers in the various areas in Figures 1 through 4 show possible area identifiers for each area.
An area map is a set of edge segments in the plane. There are three integers associated with each edge segment in the area map: a group identifier, a left area identifier, and a right area identifier. An area map is represented by a linked list of information that is constructed by calls to routines in the Areas utility. In the Fortran code, the area map (called MAP in the examples) is stored in a large integer array that is passed from subroutine to subroutine.
There are two additional important concepts in Areas. A mask is a set of areas that allows you to do one thing in one set of areas and another in the masked areas. Masks are frequently used to hide latitude and longitude lines, protect areas from contours, change line or fill style and color, and so on. It is also important to understand how areas can be created from more than one edge group.
An area map can be used as a mask to hide parts of your plot from line-drawing, color-filling, and other routines so that you draw on or color only portions of your plot. The upper figure shows a plot where latitude and longitude lines are masked over the land masses so that they are only drawn over the ocean.
The previous module describes how to define an area by starting at an
edge segment and tracing the leftmost edge segments until returning to
the first edge segment. In that case, we traced out the state of
Nevada as an area. In the lower figure, we overlay a contour edge
group on the edge group defining the map of the western United States.
Using the same approach as we used before, we now trace out the shaded
area, the portion of a contour level over Nevada. Notice that as we do
this, we choose the leftmost edge, regardless of which edge group it's
in. In this case, the area has two area identifiers, one for the
contour band that it's in and one for the state that it's in. As
explained in module "Ar 1.4
Areas definitions: Groups of edges," contour levels are in
edge group 3, and geographic maps are in edge group 1 by default.
Areas is one of the more complex utilities to use because it requires you to make at least three subroutine calls before a plot is drawn. It will probably also require you to write a subroutine or two to define exactly what is supposed to be done with specific areas after they are defined. For example, if you want to draw lines over a portion of your plot and color other portions of your plot, you would probably have to write two routines to do this: one to draw the lines in selected areas and one to color-fill other areas.
The functional outline below shows how Areas works.
--------------------------------------------------------------
* 1. Open GKS
2. Set parameters
* 3. Initialize area map
* 4. Add edges to area map
* 5. Preprocess area map
* 6. Do one or more of the following:
* Break a polyline into segments, each of which lies
entirely in a single area, then process each segment
individually
* Retrieve area identifiers for areas containing a given
point
* Scan an area map to retrieve each subarea, then
process each subarea individually
7. Debug plot
* 8. Call FRAME
* 9. Close GKS
--------------------------------------------------------------
* Steps needed to produce a plot with Areas
Areas allows the advanced user to adjust a few internal variables via its parameter setting and retrieving routines. All Areas parameters have integer values. These parameters are useful mostly for debugging Areas problems. As such, most users will not need the ARGETI and ARSETI routines, or the parameters that they address.
The first argument in a call to ARGETI or ARSETI is the parameter name. The second argument is either a variable in which the value of the named parameter is to be returned (via ARGETI) or an expression specifying the desired new value of the parameter (via ARSETI).
CALL ARGETI (PNAM, IVAL)
CALL ARSETI (PNAM, IVAL)
- ARGETI
- Retrieves integer parameter information.
- ARSETI
- Sets integer parameter information.
- PNAM
- Character expression, Input or Output---The name of the parameter that you want to set or retrieve. Two-character parameter names appear in examples and discussions in place of PNAM. The parameter names and their meanings are enclosed in single quotation marks in the Fortran code to ensure that they are treated as character expressions rather than real or integer variables. Only the first two characters within the quotation marks are examined.
- IVAL
- Integer, Input or Output---An integer value for ARSETI, or an integer variable to hold the parameter value for ARGETI.
Because most users will not need to adjust Areas parameters, this
chapter discusses only the most frequently used parameter, DB. For an
example of how to set an Areas parameter, please see module "Ar 4.1 Debugging options:
Setting parameters."
Module "Ar 1.2 Table of Areas
parameters" is a quick reference guide to all Areas
parameters. Section "Ar 5. Areas parameter
descriptions" provides more thorough descriptions of the Areas
parameters.
Areas is one of several NCAR Graphics utilities that requires initialization before use. This section covers initializing an area map, adding edges to the area map, and preprocessing the area map.
Before you can use Areas, you must initialize an area map and add edges to it. Preprocessing area maps can then be done either directly or implicitly in the routines that produce results.
The last part of this section discusses preprocessing for users who want a deeper understanding of how Areas works.
--------------------------------------------------------------
1. Open GKS
2. Set parameters
* 3. Initialize area map
* 4. Add edges to area map
* 5. Preprocess area map
6. Do one or more of the following:
Break a polyline into segments, each of which lies
entirely in a single area, then process each segment
individually
Retrieve area identifiers for areas containing a given
point
Scan an area map to retrieve each subarea, then
process each subarea individually
7. Debug plot
8. Call FRAME
9. Close GKS
--------------------------------------------------------------
* Steps discussed in this section.
Before any other Areas subroutine can be used, Areas must be initialized using ARINAM.
1 CALL COLOR
2 CALL ARINAM (MAP, LMAP)
3 CALL AREDAM (MAP, XGEO, YGEO, NMAP, 1, 2, 1)
4 CALL AREDAM (MAP, XCNTR, YCNTR1, NPTS, 3, 2, 1)
5 CALL AREDAM (MAP, XCNTR, YCNTR2, NPTS, 3, 3, 2)
6 CALL AREDAM (MAP, XCNTR, YCNTR3, NPTS, 3, 4, 3)
7 CALL AREDAM (MAP, XCNTR, YCNTR4, NPTS, 3, 5, 4)
8 CALL AREDAM (MAP, XCNTR, YCNTR5, NPTS, 3, 6, 5)
CALL ARINAM (MAP, LMAP)
- MAP(LMAP)
- Integer array, Output---An integer array in which an area map is to be constructed. Each vertex for an edge segment in your area map requires ten words in the array MAP. Remember that the total number of vertexes includes those added at each intersection of edges, and those added when long edge segments are broken into shorter edge segments.
- LMAP
- Integer, Input---Length of the MAP array. As part of the initialization process, the value of LMAP is stored in MAP(1); that way, it need not be given as an argument in calls to other routines that use the area map.
Line 1 of the carmap.f code segment sets up a color table. Line 2 initializes a very large area map. In lines 3 through 8, edges are added to the area map using AREDAM, which is discussed in the next module. In the plot, notice that each area has 4 values printed in it. These values represent the group and area identifiers for each area. The next module describes how they are set.
If you set the area map too small, the error message:
ERROR 5 IN AREDAM - AREA-MAP ARRAY OVERFLOW
occurs when edges are added to the area map. There is no good way to predict exactly how large the area map should be before adding edges to it. Generally, if you are starting with a complex geographic map, start with an area map of 150,000 words. Most contour plots need between 50,000 and 100,000 words. Once you have your area map large enough, you can determine how much space was used by calculating:
MAP(1) - MAP(6) + MAP(5)
where MAP is your area map after all the edges have been added to it. Areas routines ARPRAM and ARSCAM use more space during execution than after execution, so be sure to calculate your area map size after your last call to an Areas routine, and make your area map slightly larger than the calculated size.
- Add code to carmap.f to determine how much of the area map is used, then reduce the size of the area map to some reasonable value.
Each edge is added to the area map with a call to AREDAM. When calling AREDAM directly, you must also assign group and area identifiers for each edge.
1 CALL COLOR
2 CALL ARINAM (MAP, MAPSIZ)
3 CALL AREDAM (MAP, XGEO, YGEO, NMAP, 1, 2, 1)
4 CALL AREDAM (MAP, XCNTR, YCNTR1, NPTS, 3, 2, 1)
5 CALL AREDAM (MAP, XCNTR, YCNTR2, NPTS, 3, 3, 2)
6 CALL AREDAM (MAP, XCNTR, YCNTR3, NPTS, 3, 4, 3)
7
8 CALL AREDAM (MAP, XCNTR, YCNTR5, NPTS, 3, 6, 5)
CALL AREDAM (MAP, XCURVE, YCURVE, LCURVE, IGRP, IDLEFT, IDRGHT)
- MAP(LMAP)
- Integer array, Workspace---An array containing an area map initialized by a call to ARINAM. As part of initializing the area map, ARINAM stores the dimension of MAP in MAP(1); therefore, the dimension does not have to be given as an argument in calls to AREDAM.
- XCURVE(LCURVE)
- Real array, Input---The X coordinates, in the user coordinate system, of the points defining an edge.
- YCURVE(LCURVE)
- Real array, Input---The Y coordinates, in the user coordinate system, of the points defining an edge.
- LCURVE
- Integer, Input---The number of coordinates in XCURVE and YCURVE.
- IGRP
- Integer, Input---Identifier of the group to which the edge belongs.
- IDLEFT
- Integer, Input---Identifier of the area that lies to the left of the edge.
- IDRGHT
- Integer, Input---Identifier of the area that lies to the right of the edge.
Line 1 of the carmap.f code segment sets up a color table. Line 2 initializes the area map, and lines 3 through 8 add edges to the area map. Line 3 adds the continental outline in group 1 just as Ezmap does by default. Since the outline is drawn counterclockwise, the interior of the continent is to the left of each edge segment. The area identifier to the left of each edge segment is set to be 2, giving the interior of the continent an area identifier of 2. Similarly, the area to the right of each edge segment is given an area identifier of 1, giving the oceans area identifiers of 1.
Line 4 adds the bottom contour line. It is drawn from left to right, and it is placed in group 3 in the area map just as Conpack does by default. The area to the right of each edge segment in this contour line is given an area identifier of 1, giving the bottom contour level an area identifier of 1. The area to the left of each edge segment in the first contour line is given area identifier 2. When the second contour line is added to the area map in line 5, it is also drawn from left to right, and it is added to group 3 in the area map. The area to the right of each edge segment in this contour line is also given area identifier 2. This gives the second contour band the same area identifier for both contours so that when the area map is reconciled, the second contour level will be given area identifier 2. Lines 6 through 8 add the rest of the contour lines in a similar fashion, giving each contour band a different area identifier.
Conpack uses much the same principle as this example. Each contour line in a contour plot is put in the same group. Then area identifiers for the area above and the area below the contour are assigned, given that Conpack always draws a contour so that the area below the contour line is to the left, and the area above the contour line is to the right. Each contour level is given a single area identifier, by default starting with 1 and going to n, where n is the number of contour levels in the plot.
If you count vertexes in each edge group, you can see that the number of vertexes nearly doubles when you put both edge groups in the same area map. By putting both vertexes in the same area map, you gain control over each sub-area. However, you also increase your area map size, and you increase the CPU time needed to process each area. Since Areas is one of the most CPU-intensive utilities in NCAR Graphics, this can significantly increase the time needed to run your program. In many cases, you can produce results significantly faster by putting your geographic outlines in one area map, your contour lines in another area map, then using masking.
- In the graphic, note that each area has two area identifiers and two group identifiers. Identifiers marked in red are for geographic information, and those in blue are for contour information. Be sure you understand how each area identifier is set.
- Put each contour in a different edge group, and explain what happens.
This group of four modules can deepen your grasp of Areas, but it may be skipped because area map processing is done automatically by the Areas action routines. You call ARPRAM to give information about your edge groups to Areas; this can significantly shorten processing time. If you don't call ARPRAM, other routines will call it for you. This module tells how to call ARPRAM and briefly describes the seven processing steps it does. The next three modules elaborate on these processing steps.
Frame 1 Frame 2
Frame 3 Frame 4
CALL ARPRAM (MAP, IF1, IF2, IF3)
- MAP(LMAP)
- Integer array, Workspace---An array containing an area map that has been initialized by a call to ARINAM and to which edges have been added by calls to AREDAM.
- Note: As part of initializing the area map, ARINAM stores the dimension of MAP in MAP(1); therefore, the dimension does not have to be given as an argument in calls to ARPRAM.
- IF1
- Integer, Input---If you set IF1 nonzero, ARPRAM examines a pair of edge segments for intersection only if one of the pair has a left or right area identifier that is zero or negative. This would be appropriate for contour lines that intersect the perimeter, but that are known not to intersect each other.
- IF2
- Integer, Input---If you set IF2 nonzero, ARPRAM does not check for dangling edges. This would be appropriate for contour lines that are known not to have any such edges. This is not appropriate for Ezmap boundary lines, because the Ezmap dataset contains some edge segments (small islands) that are formed by unclosed curves.
- IF3
- Integer, Input---If you set IF3 nonzero, ARPRAM speeds up the process of adjusting area identifiers by omitting the consideration of "holes" in the areas examined. This is appropriate for contour lines as long as the left and right area identifiers provided at each contour level are consistent with those provided at the adjacent levels.
Preprocessing an area map is a seven-step process that removes "dangling" edge segments (those that do not contribute to enclosing any area) and makes all of the area-identifier information in the area map consistent:
Step 1: ARPRAM shortens edge segments whose projections on the X axis are more than twice as long as the average. ARPRAM does this by interpolating points along their lengths. This improves efficiency when executing other parts of the algorithm.
Step 2: ARPRAM locates all intersections of edge segments and interpolates the intersection points along these edge segments. This step can take a lot of time. If you set IF1<>0, then pairs of edge segments are examined for intersections only if one of the pair has a left or a right area identifier that is zero or negative.
Step 3: ARPRAM removes coincident edge segments (those with identical endpoints). If two coincident edge segments belong to the same group, one of them is removed; if they belong to different groups, one of them is modified in such a way that it will appear to be present when looking for areas defined by edge segments in a particular group, but absent when looking for areas defined by all edge segments.
Step 4: ARPRAM searches for and removes "dangling" edges (those that do not contribute to enclosing any area). This step is skipped if IF2<>0.
Step 5: ARPRAM looks for holes in the areas defined by each edge group and draws connecting edge segments between the edge and the hole. This step is skipped if IF3<>0.
Step 6: ARPRAM adjusts area identifier information in the area map. It examines all the edge segments of each area in each group to see what area identifier should be assigned to the area, and then makes adjustments.
Step 7: Connecting lines that were inserted in step 5 are removed from the area map.
You can put edges in an area map, preprocess it, add more edges, preprocess it again, and so on.
No discussion is given about the code that generates plots for the four modules in this preprocess area map group. Rather than analyzing the code, study the plots to learn how ARPRAM operates. Frames 1 through 3 from this code show the three edge groups that are added to the area map before processing, and frame 4 shows the entire area map before processing. The arrows show the direction of each edge segment, and the small numbers on each side of each edge segment are the left and right area identifiers for the segment.
The first two steps that ARPRAM does break up edge segments into shorter lengths. Step 1 shortens long edges for faster processing, and step 2 finds all points where edges intersect and breaks edge segments at those points.
Frame 1 Frame 5
Frame 8 Frame 12
In step 1, any edge segment whose projection on the X axis is longer than twice the average is replaced by two or more edge segments whose projections on the X axis are about the same as the average. Examine the edges in group 1 in frame 1 and note that the perimeter is drawn with unbroken top and bottom edge segments. In frame 5, however, you can see that those edge segments are divided in the middle by step 1. Step 1 is done on all edge groups in the area map to speed up later steps in processing the area map.
Step 2 breaks up edge segments at points of intersection with other edge segments. Compare frames 8 and 12; notice how Areas adds a vertex wherever an edge segment meets another edge segment, regardless of which groups the edge segments are in. It also shows that both edge segments are divided into shorter edge segments at that vertex.
More precisely, ARPRAM decides where to break edge segments if an edge segment (P,Q) is crossed by another edge segment at point I, where I is distinct from P and Q. ARPRAM then breaks (P,Q) into two edge segments, (P,I) and (I,Q). Clearly, this process can considerably increase the size of the area map. Also note that ARPRAM finds these points of intersection regardless of which group or groups the affected edge segments are in. Frame 12 shows the area map after this step is completed.
As discussed in the previous module, step 2 can be abbreviated by calling ARPRAM with IF1<>0. In this case, edge segments are examined for intersection only if one of the pair has a right or a left area identifier that is zero or negative. This is appropriate for contour lines, or if you know that you have equally "nice" edges that don't cross anywhere except at the perimeter.
In its next two steps, ARPRAM removes coincident (overlapping) edges and dangling edges that don't contribute areas to the area map. Step 3 removes coincident edges, and step 4 removes dangling edges.
Frame 10 Frame 14
Frame 11 Frame 19
Step 3 searches the area map for coincident edge segments (edge segments with identical end points). If both edge segments are in the same edge group, then one of the edges is removed from the area map, and area identifier information is discarded (the area identifier is set to zero) unless it is self-consistent. For example, if you put edge segment (P,Q) into the area map with a right area identifier of 1, and later give (P',Q') with the same end points a right area identifier of 2, then the area map will contain a single edge segment (P,Q) with right area identifier 0. If the coincident edge segments are in different groups, then both will appear to be present when looking for areas defined by edges in the same group, and only one will appear to be present when looking for areas defined by all edges. Examine the third diagonal line from the upper left in frames 10 and 14. In frame 10, the overlapping left and right area identifiers show the coincident edges, and in frame 14, the area identifiers for the affected segments are set to zero.
This step ensures that areas are correctly traced when tracing an area boundary using the method of following rightmost or leftmost edge segments.
Step 4 removes dangling edge segments (edge segments that do not contribute to forming any areas). A dangling edge is any edge that has an endpoint with only one edge segment in an edge group emanating from it. Examples of dangling edges include a single line segment with endpoints not on the perimeter, short curves drawn by Ezmap to represent islands, and in frame 11, the tail on the square. Note that the tail is drawn as part of group 5; however, after step 4 is done for frame 19, the tail has been removed from the area map.
The last three steps that ARPRAM does are all used to reconcile area identifier information. ARPRAM adds temporary edge segments to the area map to remove any holes, then it reconciles area identifier information and removes the temporary edge segments.
Frame 17 Frame 21
Simplified figure Frame 25
The concept of a "hole" is hard to define precisely. Concentric circles are examples of areas with holes. Holes are shown in frame 17 on the preceding page: the square forms a hole in the area whose outer boundary is the octagon, and the octagon forms a hole in the area whose outer boundary is the perimeter. Closed contour lines are also examples of holes.
Holes are important because an area with a hole does not have a boundary that is a simple closed curve. Instead, the boundary consists of two or more simple closed curves; one of these forms the outer boundary of the area and the others form the holes. This complicates things for Areas because, in attempting to reconcile contradictory area identifier information, it can only look at one closed curve at a time; contradictory information on two different parts of the boundary of an area with holes would therefore not be reconciled.
In step 5, ARPRAM identifies holes and adds temporary edge segments to the area map to, in effect, remove the holes. These edge segments ensure that all the areas in the area map have boundaries that are simple closed curves. Thus, as ARPRAM traces the boundary of one area, it sees all of the area identifiers that have been provided by the user for that area.
Frames 17 and 21 show edge group 1 before and after step 5. Consider the area having the perimeter as an outer boundary and the octagon as a hole. Information provided with the perimeter says that its area identifier should be 0, but information provided with the octagon says that its area identifier should be 1. Step 5 inserts a vertical connecting line between the perimeter and the octagon, making it possible to trace the boundary of the area as a single simple closed curve, so that all area-identifier information can be considered at once and the contradiction can be seen and reconciled.
ARPRAM methodically traces all the area boundaries defined by an area map: it repeatedly selects an edge segment to start with and makes either leftmost turns only or rightmost turns only until it arrives back at the starting point. Edge segments seen while doing this are marked so that a given boundary is only traced once.
Holes are detected by keeping track of the total change in direction along each boundary as it is traced. Net changes in angular direction to the right (clockwise) are considered negative, and net changes in angular direction to the left (counterclockwise) are considered positive. The net change along a boundary will always be exactly +360 or -360 degrees. Either of the following two cases is sufficient to confirm the existence of a hole:
- If we always take rightmost turns and yet we trace a boundary in the counterclockwise direction (the total change in angular direction is +360 degrees), then what is to the left of our path is a hole in the area to the right.
- If we always take leftmost turns and yet we trace a boundary in the clockwise direction (the total change in angular direction is 360 degrees), then what is to the right of our path is a hole in the area to the left.
An example of a hole is given in the "Simplified figure" on the preceding page. If we start at S, and we proceed in the direction of the arrowhead, and we always take rightmost turns, the total change in direction is +360 degrees (counterclockwise). Thus, the area to the left is a hole in the area to the right.
If you set IF3 nonzero, ARPRAM will not look for holes. This should only be done when you can be absolutely sure that area identifier information for holes in an area does not contradict that for the outer boundary of the area.
In step 6, ARPRAM again traces all of the boundaries in the area map. This time, ARPRAM assumes that each such boundary defines precisely one area, and it reconciles the area identifier information given for that area with the edge segments defining that boundary. After ARPRAM has picked an area identifier for the area, it replaces all the area identifiers for those edge segments with the value selected. If any of the area identifiers examined is negative, a -1 will be used; if all of the area identifiers examined are zeroes, a 0 will be used; otherwise, the nonzero value most recently entered in the area map will be used.
In step 7, ARPRAM removes the temporary edge segments that were added in step 5.
In frame 25, note that the temporary edge segments added in step 5 have been removed and that the area identifiers on the inside of the perimeter now match the values on the outside of the octagon.
After Areas is initialized, you can access and process each area in your area map individually. This flexibility allows you to do much more than the masking, filling, and line drawing discussed here. This section covers retrieving area ids and group ids, masking and drawing lines through an area map, writing a masking or line-drawing routine, filling areas, writing an area-fill routine, and describes how Areas interacts with contours and maps.
--------------------------------------------------------------
1. Open GKS
2. Set parameters
3. Initialize area map
4. Add edges to area map
5. Preprocess area map
* 6. Do one or more of the following:
* Break a polyline into segments, each of which lies
entirely in a single area, then process each segment
individually
* Retrieve area identifiers for areas containing a given
point
* Scan an area map to retrieve each subarea, then
process each subarea individually
7. Debug plot
8. Call FRAME
9. Close GKS
--------------------------------------------------------------
* Steps discussed in this section.
In the preceding modules, you may have wondered how we got the information to print out the group and area identifier arrays for each area. ARGTAI takes a point in the user coordinate system and returns information about the area containing that point. This information can then be printed with a WRITE statement or plotted on your area map using Plotchar.
1 DO 20, I=1, 12
2 CALL ARGTAI (MAP, X(I), Y(I), IAREA, IGRP, NGRPS, NAI, 0)
3 DO 30, J=1, 2
4 WRITE (STRING, 11) J, IAREA(J), J, IGRP(J)
5 IF (IGRP(J) .EQ. 1) THEN
6 CALL GSPLCI(1)
7 CALL PLCHHQ (X(I), Y(I), STRING, .01, 0., 0.)
8 ENDIF
9 IF (IGRP(J) .EQ. 3) THEN
10 CALL GSPLCI(6)
11 CALL PLCHHQ (X(I), Y(I) -.018, STRING, .01, 0., 0.)
12 ENDIF
13 30 CONTINUE
14 20 CONTINUE
CALL ARGTAI (MAP, XCD, YCD, IAREA, IGRP, ISIZ, NGRPS, ICF)
- MAP(LMAP)
- Integer array, Workspace---An array containing an area map that has been initialized by a call to ARINAM and to which edges have been added by calls to AREDAM. If you did not preprocess the area map by calling ARPRAM, ARGTAI calls it before doing anything else.
- Note: As part of initializing the area map, ARINAM stores the dimension of MAP in MAP(1); therefore, the dimension does not have to be given as an argument in calls to ARGTAI.
- XCD
- Real, Input---The X coordinate, in the current user coordinate system, of a point about which you want to obtain information.
- YCD
- Real, Input---The Y coordinate, in the current user coordinate system, of a point about which you want to obtain information.
- IAREA(ISIZ)
- Integer array, Output---The array in which area identifiers for the area containing the specified point will be returned.
- IGRP(ISIZ)
- Integer array, Output---The array in which group identifiers for the area containing the specified point will be returned.
- ISIZ
- Integer, Input---Dimension of the arrays IAREA and IGRP. ISIZ>=NGRPS.
- NGRPS
- Integer, Output---The number of values returned in IAREA and IGRP. NGRPS equals the number of groups of edges that you put in the area map.
- ICF
- Integer, Input---Flag set nonzero by you to indicate that the definition of the user coordinate system has been changed since the last call to ARGTAI; in this case, calls to GETSET must be executed by ARGTAI. If you set the flag to zero, Areas assumes that the information retrieved previously is still correct and skips calls to GETSET. This option is valuable if you are planning to make thousands of calls to ARGTAI with the same area map: do the first one with ICF=1 and the rest with ICF=0.
In the carmap example, the arrays X and Y contain the X and Y coordinates of points in the area map that are chosen to represent many of the areas. Line 1 of the carmap.f code segment sets up a loop so that we can find out area and group identifier information for each point. Line 2 retrieves the area and group identifier information. Since we know that geographic area identifiers are in group 1, lines 5 through 8 write them out in red, and lines 9 through 12 write the contour area identifier information out in blue.
ARGTAI has also been used to determine whether a specific point is over land or water. A user created a fine grid over the globe and used ARGTAI to determine if each point in the grid was over land or water. Then he used the resulting array to tell his model which points were over land and which points were over water.
ARDRLN takes a polyline (or curve) and an area map and breaks the polyline into pieces, each of which lies entirely within a single area defined by the area map. It calls user routine LPR once for each piece of the polyline. You can then program LPR to do whatever you want with the pieces. For example, you can change line color, dash pattern, or mask portions of the line that lie over selected areas. ARDRLN is most frequently used in conjunction with Ezmap to draw latitude and longitude lines over the oceans.
1 EXTERNAL MASK
2 CALL ARINAM (MAP, LMAP)
3 CALL AREDAM (MAP, XGEO, YGEO, NGEO, 1, 2, 1)
4 CALL ARDRLN (MAP, XCNTR, YCNTR1, NPTS, XWRK, YWRK, NWRK,
+ IAREA, IGRP, NGRPS, MASK)
5 CALL ARDRLN (MAP, XCNTR, YCNTR2, NPTS, XWRK, YWRK, NWRK,
+ IAREA, IGRP, NGRPS, MASK)
6 CALL ARDRLN (MAP, XCNTR, YCNTR3, NPTS, XWRK, YWRK, NWRK,
+ IAREA, IGRP, NGRPS, MASK)
7 CALL ARDRLN (MAP, XCNTR, YCNTR4, NPTS, XWRK, YWRK, NWRK,
+ IAREA, IGRP, NGRPS, MASK)
8 CALL ARDRLN (MAP, XCNTR, YCNTR5, NPTS, XWRK, YWRK, NWRK,
+ IAREA, IGRP, NGRPS, MASK)
CALL ARDRLN (MAP, XCD, YCD, NCD, XCS, YCS, MCS, IAREA, IGRP,
+ ISIZ, LPR)
- MAP(LMAP)
- Integer array, Workspace---The area-map array.
- XCD(NCD)
- Real array, Input---The X coordinates, in the user coordinate system, of points defining the curve.
- YCD(NCD)
- Real array, Input---The Y coordinates, in the user coordinate system, of points defining the curve.
- NCD
- Integer, Input---Number of points for the curve.
- XCS(MCS)
- Real array, Workspace---Used in calls to LPR to hold the X coordinates of the curve segment contained in a given area.
- YCS(MCS)
- Real array, Workspace---Used in calls to LPR to hold the Y coordinates of the curve segment contained in a given area.
- MCS
- Integer, Input---Dimension of each of the arrays XCS and YCS.
- IAREA(ISIZ)
- Integer array, Workspace---The area identifier array used in calls to LPR.
- IGRP(ISIZ)
- Integer array, Workspace---The group identifier array used in calls to LPR.
- ISIZ
- Integer, Input---Dimension of the each of the arrays IAREA and IGRP. ISIZ>=n, where n is the number of groups in the area map.
- LPR
- Subroutine---A line-processing routine that you supply for ARDRLN to call. LPR must be declared external in the routine that calls ARDRLN. See the next module.
ARDRLN overlays the curve that you give it onto the area map. It breaks the curve into pieces, one for each area that it crosses. LPR is then called with each piece of the curve. For the green contour in the carline example, ARDRLN breaks the contour into five pieces. The first piece is from the left edge to Mexico. ARDRLN calls the line-processing routine MASK, which checks to see if this piece is over water; if it is, ARDRLN draws it. ARDRLN then calls MASK with the second piece (over Mexico) and does not draw it because it is not over water. ARDRLN calls MASK with the third piece over the Gulf of Mexico and draws it. It does not draw the fourth piece over Central and South America, then it draws the fifth piece over the Atlantic Ocean. This shows how you can make a line-processing routine do something different with each piece of a curve.
Before executing the first call to LPR, ARDRLN calls GETSET to retrieve the current user-system mapping parameters, then executes the statement:
CALL SET (VPL, VPR, VPB, VPT, VPL, VPR, VPB, VPT, 1)
where VPL, VPR, VPB, and VPT are the viewport left, right, bottom, and top coordinates in NDCs.
This ensures correct results if the NDCs in XCS and YCS are used in calls to such routines as GPL and CURVE, and this allows clipping at the edges of the viewport. LPR may make its own SET call to achieve some other effect. Before returning control to the calling routine, ARDRLN calls SET again to restore the original mapping parameters.
In the carline example, XGEO and YGEO are arrays containing the X and Y coordinate information for the "geographical" outline. XCNTR and YCNTR1-4 are the X and Y coordinate arrays for the "contour" lines. Line 1 declares our masking routine to be external so that Fortran doesn't interpret the subroutine name as a variable. Line 2 initializes Areas. Line 3 adds the geographic outline to the area map in group 1 for consistency with Ezmap defaults. Lines 4 through 8 show examples of calling ARDRLN to mask the contour lines over land masses, where MASK is the actual line-drawing and masking routine.
When drawing lines with ARDRLN, it is best not to add the lines to the area map.
By itself, Areas won't actually do anything after it defines areas that can be acted on. You have to write one or two routines that actually do something because only you know what you want done inside the defined areas. The following subroutine is an example that checks which area the line is in and decides what color the line should be before drawing it. It also decides when not to draw the line.
1 SUBROUTINE MASK (XC, YC, MCS, IAREA, IGRP, NGRPS)
2 INTEGER IAREA(NGRPS), IGRP(NGRPS), IDGEO
3 REAL XC(MCS), YC(MCS)
4 IDGEO = -1
5 DO 10, I=1, NGRPS
6 IF (IGRP(I) .EQ. 1) IDGEO = IAREA(I)
7 10 CONTINUE
8 IF (IDGEO .EQ. 1) THEN
9 CALL CURVED (XC, YC, MCS)
10 ENDIF
11 RETURN
12 END
SUBROUTINE LPR (XCS, YCS, NCS, IAREA, IGRP, ISIZ)
DIMENSION XCS(*), YCS(*), IAREA(*), IGRP(*)
- LPR
- Subroutine---An area map line-processing routine you supply for ARDRLN to call. LPR must be declared external in the routine that calls ARDRLN.
SUBROUTINE LPR (XCS, YCS, NCS, IAREA, IGRP, ISIZ)
DIMENSION XCS(*), YCS(*), IAREA(*), IGRP(*)
- (code to process polyline defined by XCS, YCS, IAREA, and IGRP)
RETURN
END
- XCS(NCS), YCS(NCS)
- Real arrays, Input---Hold the X and Y coordinates, in NDCs, of NCS points defining a piece of the original polyline that is completely contained in one of the areas defined by the area map.
- NCS
- Integer, Input---Number of X and Y coordinates in the arrays XCS and YCS.
- IAREA(ISIZ), IGRP(ISIZ)
- Integer arrays, Input---Hold ISIZ pairs of area identifiers and group identifiers, respectively, for the area in which the piece of the polyline lies. For each value of I from 1 to ISIZ, IAREA(I) is the area identifier for the area with respect to the group of edges specified by the group identifier IGRP(I).
- ISIZ
- Integer, Input---Number of values given in IAREA and IGRP. ISIZ>=n, where n is the number of groups of edges that you put in the area map.
If the group identifier is 1, then we know the following:
- The area outside the perimeter has area identifier -1
- The areas over ocean have area identifier 1
- The areas over land have area identifier 2
Since we want to draw the contour lines only over ocean in this case, we first retrieve the area identifier for group 1 in lines 5 through 7 of the carline.f code segment. Because there is more than one group, it is necessary to check each group identifier individually in a DO loop. If the area identifier for group 1 is 1, then we draw the contour line segment in line 9.
- Using the carline example, mask the areas over water and not land.
- Draw a vertical line down the middle of the plot, and give a different color to the segment passing through each contour level.
Given an area map, ARSCAM scans it from left to right, extracting the definitions of all the areas in the area map. ARSCAM then processes each area using the user routine APR. Since you supply APR, you can use ARSCAM to do almost anything with the areas. ARSCAM is most frequently used to color-fill or pattern-fill regions.
1 EXTERNAL FILL
2 CALL GSFAIS (1)
3 CALL ARINAM (MAP, LMAP)
4 CALL AREDAM (MAP, XGEO, YGEO, NMAP, 1, 2, 1)
5 CALL AREDAM (MAP, XCNTR, YCNTR1, NPTS, 3, 2, 1)
6 CALL AREDAM (MAP, XCNTR, YCNTR2, NPTS, 3, 3, 2)
7 CALL AREDAM (MAP, XCNTR, YCNTR3, NPTS, 3, 4, 3)
8 CALL AREDAM (MAP, XCNTR, YCNTR4, NPTS, 3, 5, 4)
9 CALL AREDAM (MAP, XCNTR, YCNTR5, NPTS, 3, 6, 5)
10 CALL ARSCAM (MAP, XWRK, YWRK, NWRK, IAREA, IGRP, NGRPS, FILL)
CALL ARSCAM (MAP, XCS, YCS, MCS, IAREA, IGRP, ISIZ, APR)
- MAP(LMAP)
- Integer array, Workspace---An array containing an area map that has been initialized by a call to ARINAM and to which edges have been added by calls to AREDAM. If you did not preprocess the area map by calling ARPRAM, ARSCAM calls it before doing anything else.
- Note: As part of initializing the area map, ARINAM stores the dimension of MAP in MAP(1); therefore, the dimension does not have to be given as an argument in calls to ARSCAM.
- XCS(MCS), YCS(MCS)
- Real arrays, Workspace---Hold the X and Y coordinates defining the edge of a given area, for use by ARSCAM in calls to the area-processing routine APR.
- MCS
- Integer, Input---Dimension of each of the arrays XCS and YCS.
- IAREA(ISIZ), IGRP(ISIZ)
- Integer arrays, Workspace---Hold ISIZ pairs of area identifiers and group identifiers, respectively. Used by ARSCAM in calls to the area-processing routine APR.
- ISIZ
- Integer, Input---Dimension of each of the arrays IAREA and IGRP. ISIZ>=n, where n is the number of groups in the area map.
- APR
- Subroutine---An area-processing routine that is called by
ARSCAM and that must be declared external in the routine that calls
ARSCAM. You must supply the routine APR; it is described in detail in
the next module.
When ARSCAM is called by carfill.f, it locates the contour band at the bottom of the plot located over the Pacific Ocean, and passes that area to the area-processing routine FILL. Because the area is over the ocean, FILL colors the area blue; it then returns to ARSCAM. The next area ARSCAM locates and passes to FILL is the portion of the contour band directly above the first. FILL colors the area blue and returns to ARSCAM. After ARSCAM has passed all the contour bands over the Pacific to FILL, the next area it finds is over the western part of Alaska. ARSCAM passes this area to FILL, which recognizes that the area is in the fifth contour band and over land, so it fills the area with the color specified for the fifth contour band over land: magenta. Next, ARSCAM finds the area over eastern Alaska and northwest Canada and passes it to FILL. FILL recognizes that the area is in the sixth contour band and over land, so it fills it with white. ARSCAM then locates the area over central Canada and passes it to FILL. FILL recognizes that the area is in the fifth contour band and over land, so it fills it with magenta, and so on.
All regions having the same group and area identifiers are treated the same way. In this example, every separate piece of a single contour band over land is filled with the same color. This is why the color of western Alaska matches the color of central Canada. If you want to color a contour over Canada differently than the same contour over the US, you would have to add the political boundaries to your area map, then give Canada a different area identifier than the US. Fortunately, Ezmap can do this for you.
Line 1 of the carfill.f code segment declares our filling
routine FILL to be external so that when we call it, Fortran doesn't
interpret the subroutine name as a variable and then core dump. Line 2
sets the GKS fill style to be solid because the default is
"hollow" fill. Line 3 initializes Areas. Lines 4 through 9
add "geographic" and "contour" lines to the area
map using group 1 for geographic outlines and group 3 for contour
lines; this ensures consistency with Conpack and Ezmap defaults. Line
10 does the area fill by calling ARSCAM, which examines the area map
and feeds the appropriate information to the fill routine FILL. The
next module describes the FILL routine.
The Areas utility defines areas that can be acted on. You must decide what to do with these areas and then supply a routine that takes the appropriate action on the area or the line going through an area. Below is an example subroutine that supplies color and fill type for each area defined by an area map.
1 SUBROUTINE FILL (XWRK, YWRK, NWRK, IAREA, IGRP, NGRPS)
2 INTEGER IAREA(NGRPS), IGRP(NGRPS), IDGEO, IDCNTR
3 REAL XWRK(NWRK), YWRK(NWRK)
4 IDGEO = -1
5 IDCNTR = -1
6 DO 10, I=1, NGRPS
7 IF (IGRP(I) .EQ. 1) IDGEO = IAREA(I)
8 IF (IGRP(I) .EQ. 3) IDCNTR = IAREA(I)
9 10 CONTINUE
10 IF (IDGEO .EQ. 1) THEN
11 CALL GSFACI(7)
12 CALL GFA (NWRK, XWRK, YWRK)
13 ELSE IF ((IDGEO .EQ. 2) .AND. (IDCNTR .GE. 1)) THEN
14 CALL GSFACI(IDCNTR)
15 CALL GFA (NWRK, XWRK, YWRK)
16 ENDIF
17 RETURN
18 END
SUBROUTINE APR (XCS, YCS, NCS, IAREA, IGRP, ISIZ)
DIMENSION XCS(*), YCS(*), IAREA(*), IGRP(*)
- APR
- Subroutine---An area-processing routine that is called by ARSCAM and must be declared external in the routine that calls ARSCAM. You must supply the routine APR, and it must have the following form:
SUBROUTINE APR (XCS, YCS, NCS, IAREA, IGRP, ISIZ)
DIMENSION XCS(*), YCS(*), IAREA(*), IGRP(*)
- (code to process the area defined by XCS, YCS, IAREA, and IGRP)
RETURN
END
- XCS(NCS), YCS(NCS)
- Real arrays, Input---The X and Y coordinates, in NDCs, of NCS points defining a polygonal area. The last of these points is a duplicate of the first. Holes in an area are traced in a way that maximizes the probability of hardware fill working properly.
- NCS
- Integer, Input---Number of X and Y coordinates in the arrays XCS and YCS.
- IAREA(ISIZ), IGRP(ISIZ)
- Integer arrays, Input---Hold ISIZ pairs of identifiers for the area defined by XCS and YCS. For each value of I from 1 to ISIZ, IAREA(I) is the area identifier for the area with respect to the group of edges specified by the group identifier IGRP(I).
- ISIZ
- Integer, Input---Number of values returned in IAREA and IGRP. ISIZ>=n, where n is the number of groups of edges that you put in the area map.
In the carfill example, we need to distinguish between both geographic areas and contour levels since we'll only be filling contour levels over land. For this reason, lines 4 and 5 of the carfill.f code segment set our area identifier storage variables to a negative number; this avoids using information from a previous call to FILL. Lines 6 through 9 check every value in the group identifier array; if the group identifier is 1, then we store the associated area identifier in IDGEO. If the group identifier is 3, we store the associated area identifier in IDCNTR. Because of the way that Areas works, there is no way to guarantee that the first value in the area identifier array will be from your first group, the second value in the array will be from your second group, and so on.
After the area identifiers have been obtained, line 10 checks to see if we're over the ocean. If we are, then lines 11 and 12 fill the area in aquamarine. Line 13 checks to see if we're over land, and if so, then lines 14 and 15 fill the contour band with a color based on its area identifier.
APR can be used to do almost anything you want with each area in the area map. However, it is rarely used for anything except filling areas. If you are unable to produce color hardcopy output, you may want to write your fill routine so that it calls the Softfill utility.
Before executing the first call to APR, ARSCAM executes the statement:
CALL SET (VPL, VPR, VPB, VPT, VPL, VPR, VPB, VPT, 1)
where VPL, VPR, VPB, and VPT are the viewport left, right, bottom, and top coordinates in NDCs. This ensures correct results if the NDCs in XCS and YCS are used in calls to such routines as GFA and FILL, and this allows clipping at the edges of the viewport. APR may make its own SET call to achieve some other effect. Before returning control to the calling routine, ARSCAM calls SET again to restore the original mapping parameters.
- Fill the contour bands over water in different colors, and fill the land areas in white.
Because Areas is used with various other packages, and because you may not always have a complete understanding of how edge groups and area identifiers are assigned in these different cases, it is often helpful to debug a problem by drawing edge segments and area identifiers for one or more groups of edges. This section describes Areas debugging tools that can help in these situations.
--------------------------------------------------------------
1. Open GKS
* 2. Set parameters
3. Initialize area map
4. Add edges to area map
5. Preprocess area map
6. Do one or more of the following:
Break a polyline into segments, each of which lies
entirely in a single area, then process each segment
individually
Retrieve area identifiers for areas containing a given
point
Scan an area map to retrieve each subarea, then
process each subarea individually
* 7. Debug plot
8. Call FRAME
9. Close GKS
--------------------------------------------------------------
* Steps discussed in this section.
Areas currently has one debugging parameter, the debug plots flag DB. DB turns on debugging, and this directs ARPRAM to produce a series of plots that show all the edges belonging to a specified group in the area map at selected break points.
1 CALL ARINAM (MAP, MAPSIZ)
2 CALL ARSETI ('DB - DEBUG PLOTS FLAG', 1)
3 CALL AREDAM (MAP, X1, Y1, NPTS, 1, 1, 0)
4 CALL AREDAM (MAP, X2, Y2, NPTS, 1, 2, 1)
5 CALL AREDAM (MAP, X3, Y3, NVERT, 1, 3, 2)
6 CALL ARSCAM (MAP, XC, YC, MCS, AREAID, GRPID, IDSIZE, FILL)
CALL ARSETI ('DB', idb)
- DB
- Integer---The DeBug plots flag is available to help you understand an area map at various stages of the Areas run.
- 0 Debugging is off. This is the default.
- >=1 Debugging is on. At selected break points, ARPRAM produces plots showing all edge segments in the area map that have the group identifier specified by DB.
Line 1 of the cardb1.f code segment initializes Areas. Line 2 turns on debugging. ARDBPA is then automatically called to produce a plot at each of seven different break points:
- The start of ARPRAM
- After breaking up long edge segments
- After finding points of intersection
- After removing coincident segments
- After removing unclosed edges
- After looking for holes
- After updating area identifiers
Lines 3 through 5 add edges to our area map. Line 6 forces a call to ARPRAM, during which the calls to ARDBPA occur, generating the seven frames listed above. This occurs before ARSCAM fills the areas, generating an eighth frame, which is the final plot.
In each of the debug plots, notice the dots on each side of each line segment. These are numbers that represent the left and right area identifiers for that edge. To view them, you need to use the zoom feature in idt.
The area identifiers are drawn very small so that each edge segment in an area map can be marked with both left and right area identifiers.
- Turn debugging on in the carline example, and see if you understand each step in the process.
The routine ARDBPA produces a plot showing all the edges in an area map with a specified group identifier. In complex area maps, this allows you to debug a problem, group by group.
1 CALL ARINAM (MAP, LMAP)
2 CALL AREDAM (MAP, XGEO, YGEO, NMAP, 1, 2, 1)
3 CALL AREDAM (MAP, XCNTR, YCNTR1, NPTS, 3, 2, 1)
4 CALL AREDAM (MAP, XCNTR, YCNTR2, NPTS, 3, 3, 2)
5 CALL AREDAM (MAP, XCNTR, YCNTR3, NPTS, 3, 4, 3)
6 CALL AREDAM (MAP, XCNTR, YCNTR4, NPTS, 3, 5, 4)
7 CALL AREDAM (MAP, XCNTR, YCNTR5, NPTS, 3, 6, 5)
8 CALL ARDBPA (MAP, 3, 'Crossing Contours')
CALL ARDBPA (MAP, IGRP, LABEL)
- MAP(LMAP)
- Integer array, Workspace---An array containing an area map that has been initialized by a call to ARINAM and to which edges have probably been added by calls to AREDAM.
- Note: As part of initializing the area map, ARINAM stores the dimension of MAP in MAP(1); therefore, the dimension does not have to be given as an argument in calls to ARDBPA.
- IGRP
- Integer, Input---The group identifier of the group that you want to examine.
- LABEL
- Character expression, Input---The label you want on the plot.
Line 1 of the cardb2.f code segment initializes the area map. Line 2 puts our "Americas map" into group 1. Lines 3 through 7 add the contour lines to group 3. Notice that in this case, we have chosen contour lines that cross. Line 8 calls ARDBPA to examine edge group 3, which includes the crossing contours. Use the zoom feature in idt to resolve the small dots near the arrows on the contour lines into left and right area identifiers.
In the zoomed portion of the picture, notice how the left and right area identifiers are in conflict in the areas above and below the crossed contour lines. This situation occasionally arises when using smoothing in Conpack. When Areas is used in a situation like this, the contour bands may be filled in an undesired fashion, since the areas do not have well-defined identifiers.
By default, each edge in the plot appears in one of four different colors, depending on whether the area identifiers to the left and right are less than or equal to zero or greater than zero, as follows:
--------------------------------
Color Left area Right area
identifier identifier
--------------------------------
Magenta <=0 <=0
Yellow <=0 >0
Cyan >0 <=0
White >0 >0
--------------------------------
You can adjust these colors with the DC parameter.
In some cases, you may notice gray lines in your plot. This means that the same edge occurs in more than one group. In all but one of those groups, Areas negates the group identifier for the edge in the area map. This allows Areas to include the edge when it is looking at a particular group (as in ARPRAM), but omit it when it is looking at the union of all the groups (as in ARSCAM).
If ARDBPA is used for a complicated area map, the amount of output can be very large.
- Using the cardb2 example and ARDPA, examine the area identifiers for edge group 1. Why is the Americas map drawn in cyan?
Areas currently supports the following parameters.
- AT
- Integer---Arithmetic Type flag. AT allows the user to specify the type of arithmetic to be used by Areas or to find out what type it used.
- <=0 Areas decides what sort of arithmetic to use.
- 1 Real arithmetic is to be used.
- 2 Double precision arithmetic is to be used.
- 3 Multiple precision arithmetic is used; Areas should choose the base value.
- >=4 Multiple precision arithmetic is used and AT specifies the base value. (For example, the value 100 would specify the use of base-100 multiple-precision integer arithmetic.)
- Default: 0
- Note: You should use a nonzero value of AT only on the recommendation of an NCAR consultant.
- DB
- Integer---DeBug plots flag. DB allows the user to study the area map.
- 0 Debugging turned off. This is the default.
- >=1 Debugging is on. At selected break points, ARPRAM produces plots showing all edge segments in the area map that belong to the group with group identifier DB.
- DC
- Integer---Debug Colors index. DC tells ARDBPA to use color indices DC+1 through DC+5 for debugging colors.
- By default, DC=100, so ARDBPA defines and uses color indices 101 through 106.
- DI
- Integer---Direction Indicator flag. DI tells the user the direction in which the edge of an area is given. For retrieval only, the value of DI cannot be set by the user. The two possible values of DI are:
- 1 Edge of the area is given in counterclockwise order (with the interior to the left).
- 2 Edge of the area is given in clockwise order (with the interior to the right).
- Note: These values are only meaningful when the ARGETI call originates from the user-written routine (dummy argument "APR") that is called by ARSCAM. It gives you information about the area whose polygonal boundary is defined by the values of APR's arguments.
- LC
- Integer---Largest Coordinate. LC lets you specify the largest coordinate allowed in an area map. X and Y coordinates in an area map are represented by integers in the range from 0 to LC, inclusive; the default value of LC is 1000000.
- The minimum value allowed for LC is 1000; attempting to set LC<1000 gives it the value 1000. The value of LC must not be greater than the largest integer on the machine on which Areas is running; its value must also be exactly representable as a real number on that machine.
Previous Chapter Tutorial Home Next Chapter