This chapter provides the information you will need to access and use McIDAS-X data. You'll learn:
This chapter is organized into these sections:
Reading and writing data to and from disk is fundamental to many applications programs. A McIDAS-X disk file stores information that applications can randomly access by byte address using standard system library calls.
McIDAS-X disk file utilities have several characteristics that distinguish them from other file system utilities.
The table below describes the McIDAS-X library functions that you will use when programming with disk files. In McIDAS-X, a disk file is called an LW (Large Word) array file. You will notice that many of the Fortran APIs below begin with the letters lw .
Each function is described in the sections below, along with sample code illustrating its use.
The McIDAS-X library provides the Mcread and Mcwrite functions for reading and writing disk files in C. The comparable functions in Fortran are lbi and lbo , which are shown in the code fragment below. The lwi and lwo functions read and write words (4-byte groups) from a disk file.
If your application uses Fortran or C library functions for disk I/O, you must use the volnam or Mcpathname function to convert the name of the disk file into a fully qualified, system pathname/filename . This pathname is essential for locating and working on a file.
The code fragment below illustrates the use of volnam with a Fortran OPEN statement. The maximum number of characters allowed for the fully-qualified pathname is stored in the constant MAXPATHLENGTH in the Fortran INCLUDE file, fileparm.inc . The limit for C is in the file /usr/include/sys/limits.h .
Use the function lwfile to determine if a file with a given name already exists. If the file does not exist, lwfile returns a zero; it does not create the file for you. The following code fragment illustrates how to use lwfile .
c --- find out if file 'mystuff' exists status = lwfile('mystuff') if (status .ne. 0) then call sdest('The file is there!',0) else call edest('The file is NOT there!!!',0) endif ... |
To delete the contents of a disk file but not the file itself, use the lwtrunc or Mctruncate function. These functions remove the contents of the file, leaving one word (4 bytes) of 0x80808080 at the beginning of the file.
Deleting only the contents of a file and not the file itself is important if the file name appears in more than one location in the MCPATH tree. For example, if you delete the file ABC from a directory where you have write permissions but the file also exists further down your MCPATH in a directory where you have only read permissions, you won't be able to create a writable file by that name. If you use lwtrunc or Mctruncate to delete only the file's contents, you can write into it again in the future, since it resides in the original, writable directory.
The code fragment below uses the lwtrunc function to delete the contents of a file.
Use the lwd or Mcremove functions to delete a disk file, as shown in the code fragment below.
As long as a user has read and write permissions, the McIDAS-X disk file I/O subsystem will open all files and permit simultaneous access to these files for both reading and writing. If you write applications that update information in disk files, you must synchronize access to the file to avoid potential file collisions.
McIDAS-X has a locking mechanism called lock for coordinating file updates between applications or between copies of the same application. The lock function acquires exclusive use of a unique lock. If your application is using the lock, another application trying to lock the file must wait until you free the lock with the unlock function before using it.
The lock and unlock functions do not prevent other applications from reading a file or writing to it. They simply ensure an orderly means of updating a file without losing or overwriting information.
The code fragment below shows how to lock and unlock a file to protect its integrity during updating.
Use the lwcopy function to copy one disk file into another, as shown in the code fragment below.
McIDAS-X images are typically composed of atmospheric and oceanographic data, which are measured from remote sensing platforms such as satellites and radar. Images may contain any data that can be represented in a two-dimensional matrix.
Image data has several unique attributes that determine how an image is displayed on the McIDAS-X Image Window. The resolution of the image, the size of the data points and the number of spectral bands all influence the final product.
A satellite observes features by scanning small slices of the earth's surface with each pass of the sensors. The geographic width of each image line helps determine the size of the smallest surface feature the satellite can detect. This concept is called resolution .
Image resolution refers to the number of satellite image lines represented in each data point of an image line. If the image resolution is one, the image is stored at full resolution . This means that one image data point represents one satellite sensor data point. If the image resolution is four, either sampling or averaging was performed on the image so that one data point in the image represents 16 satellite scan data points. Each satellite has its own scan resolution, so an image resolution of one will mean different geographic resolutions from one satellite to another.
When copying or displaying an image, a user can modify its resolution. Image resolution can be artificially increased, or blown up , by replicating data point values, much like enlarging a 3 x 5 photograph to 8 x 10. Image resolution can be decreased, or blown down, by sampling or averaging the image data points. For example, if a blowdown factor of two is applied to a McIDAS-X image, every other data point along the line and every other line in the image is dropped out. Each data point on the displayed image represents four data points from the original, scanned image.
Data-point size refers to the number of bytes needed to accurately represent the value of a data point. Although the size varies among sensor sources, data is one, two or four bytes. For example, Meteosat-7 data is 1-byte per data point while GOES-12 data is 2-byte.
A spectral band is the wavelength in which a scanning instrument measures data; for example, band 4 for the GOES-8 Imager senses 10.7 micron wavelength radiation. These wavelengths are specific to the measuring instrument. Most satellites can measure radiation from many wavelengths simultaneously over the same geographic location.
|
For more information about the bands for satellite imagery, see Appendix B, Satellite Information . |
Image data is quantitatively useless unless it is transformed into physical units (calibration) and oriented relative to time and physical space (navigation). In addition, it is often necessary to know when and how an image was collected and processed. The actual image, along with these ancillary data , is collectively called an image object . Each image object in McIDAS-X is composed of the following blocks of information:
The API functions and the procedures for reading, writing and deleting image objects follow.
Most applications for reading image objects will do one of the following:
The table below lists alphabetically the McIDAS-X library functions for performing these tasks.
These functions are described below along with sample code.
|
See the online man pages provided with the McIDAS-X software for detailed information about any of the API functions discussed in this section. |
In ADDE, a client has the ability to request only directory blocks from a server. This allows the client to sample information on a server without transferring large amounts of image data that may not be necessary for an application. For example, the McIDAS-X IMGLIST command reads only directory blocks.
The ADDE interface to the image directory is through mcadir , which opens a connection based on a set of selection clauses for a given dataset name. The valid selection conditions are provided in the table below.
|
Selection clause |
Description |
|---|---|
|
appends center latitude/longitude, resolution, and calibration types to the comment block (default=YES) |
|
Selection clauses can restrict the search based on the image day, image start time and SSEC sensor source number. With the exception of AUX, you must specify these selection clauses as a range of values.
Note that in a mcadir call the dataset name does not contain the position field. The position, if known, is specified with the SUBSET selection clause. For example, if the application requires the first 10 images for a given dataset, the selection clause is SUBSET 1 10. If the application requires all the image directories in a dataset, the selection clause is SUBSET ALL.
If you include the selection clause AUX YES in the condition list, the image object directory server will append comment entries describing the latitude and longitude of the center element of the image, the earth area (in latitude and longitude resolution) covered by the center element of the image, and the valid calibration types for the image. These values can subsequently be parsed out with the mcpnav and mcpcal functions, as shown in the code fragment on the next page.
Once the connection is opened with mcadir , the application makes repeated calls to mcadrd until all image object directory blocks are retrieved. The mcadir call must precede the call to mcadrd , as shown below.
The directory block contains a list of ancillary information about the image. The entries in the directory block are described in the table below.
IMPORTANT: When using mcadir to read the dataset, a 65 word block is returned, which includes the absolute position number. When using mcaget, the directory is 64 words long, including only words 1-64 in the the table below.
|
Entry |
Description |
|---|---|
|
actual image start time, hhmmss; in milliseconds for POES data |
|
|
For more information about the directory block, see the AREAnnnn data structure in Chapter 6, Format of the Data Files. |
In ADDE, a client can request an entire image object, including the directory, data, line prefix, navigation, calibration and comment blocks, from a server. For example, the McIDAS-X IMGDISP and IMGCOPY commands typically request entire image objects.
To open a connection to read an image object's data block, an application must perform these two steps:
Applications use selection clauses to specify the spatial, temporal and spectral limits of the transaction with the server. This eliminates the need for the application to scan the dataset for a particular image object and also limits the size of the transmission to only that needed. The number and format of selection clauses are strictly regulated. Below is a list of the valid selection clause formats for the mcaget interface, which passes the request for an image sector from the client to the server. Additional information for each selection clause follows.
|
Selection clause |
Description |
|---|---|
|
includes the documentation from the line prefix (default=NO) |
|
|
sets the coordinate type and the coordinate positions relative to the coordinate type (default=AU 0 0) |
|
AUX YES
-- Data server: use this clause to insert the unit and scale factor of the image data into entries 58 and 59 of the directory block.
-- Directory server: use this clause to get center point and calibration information from the server as additional comment cards.
BAND --use this clause to identify a spectral band of image data. Specify BAND ALL to return all spectral bands in the source image.
CAL (obsolete) --use this clause only when the source data type is TIRO. Use CAL QTIR for quick calibration.
DAY --use this clause to specify the day of the source. Using DAY implies that the image object position (POS) is not specified in the connection request. If POS is specified, the DAY clause is ignored.
DOC YES --use this clause to include the line prefix's documentation when reading an image object.
LOCATE --use this clause to position the request spatially. LOCATE specifies a reference point from which sector bounds are determined. To specify the reference point of the image sector, use one of three coordinate systems (I=satellite, A=array or file, E=earth) with one of two standard offsets (U=upper-left, C=center); for example, EC=earth center, IU=satellite upper-left. Following the reference point are two values (ycor and xcor ) that identify the absolute position of the reference point in the chosen coordinate system. For earth coordinates, the values are latitude and longitude; for satellite coordinates, the values are line and element; for array coordinates, the values are array row and column. If you don't specify LOCATE, the default is a satellite upper-left (IU) reference point of the first scan line and pixel of the source image object.
MAG --use this clause to specify the resolution magnification factor for the line and element dimensions of the image object. Enter a negative integer for a blowdown; enter a positive integer for a blowup. For example, a magnification factor of -4 will use every fourth data point on the line and every fourth line in the image. A magnification factor of 16 for both line and element dimensions will duplicate each data point in the image 256 times (16 x 16). The mcaget function performs all line and element duplication. The application must reduce the size of the data request to allow for a blowup factor. To use a magnification factor other than one, modify the input for SIZE accordingly. For example, to get an image that is 500 x 1000 with a MAG value of 2 2, you must specify SIZE 250 500. You cannot store images with non-integer line/element resolutions in the McIDAS-X image object data structure.
POS --use this clause to specify the position of the image in the dataset. If POS isn't specified, the server uses the most current (time-relative) image. Specify POS as either absolute or time relative. Numbers less than one imply time-relative position, with zero as the most current image.
SIZE --use this clause to specify the image line and data element limits of the transaction. Specify SIZE ALL to read the entire source image object.
SU --use this clause to specify the name of an image data stretching table. These tables are generated with the McIDAS-X SU command.
TIME --use this clause to specify a range of image start times identifying specific images in a dataset. Using this clause implies that the image object position (POS) is not part of the connection request. If a POS clause does exist, the TIME clause is ignored.
|
Server-specific selection clause |
Description |
|---|---|
|
Coverage; accesses data for specified time range (realtime POES server) |
|
Once the selection conditions for the image sector are defined, the mcaget function passes the request for the image sector from the client to the server. The return status shows if the connection is open and if the request can be completed.
The mcaget function lets the application specify the units and format of the data points. Since these parameters are necessary to any data transaction, they are specified as separate parameters to the mcaget function and are not entered as selection clauses. Units may be any unit identifier valid for the image object type. A list of valid unit identifiers is available to an application through the mcadir function by specifying the AUX YES selection clause. Use the format parameter to specify the bytes-per-data point in the return array. The valid formats are I1 (1 byte per data point), I2 (2 bytes per data point) and I4 (4 bytes per data point).
If you request that 2-byte data be returned as a 1-byte representation, without going through a calibration process that converts the 2-byte data to 1-byte, the server will truncate the least significant byte.
Because mcaget requires many selection conditions to request image data, the mcasort function can be called to translate command line keyword parameters into equivalent mcaget selection clauses. Any application-level program may call mcasort to retrieve these keywords and return them as mcaget selection clauses. The table below lists the keywords for accessing image objects and their equivalent selection clauses.
If an image sector in the dataset satisfies the client request, a connection is established between the server and the client and the transaction proceeds. Because applications manipulating image data must sometimes sample data from different sources simultaneously, ADDE allows an application to receive data from multiple datasets in any order required for that application. When a request for data is initiated with mcaget , a handle identifying the request is returned. This handle is subsequently used to retrieve each line of data for the request. If your application requires accessing data from multiple requests simultaneously, you only have to manage a handle for each request.
This section describes the functions you should use to get the requested image data from the server and to read the image object's line prefix, navigation, calibration and comment blocks.
The mcalin function reads the requested image data from the server one line at a time. When all the data is read, the connection is closed. Then the mcafree function readies the environment for the next request by freeing the image handle and the memory allocated to store the image data from the previous request.
If you request multispectral data, each data point in the line will contain the measurements for the bands arranged consecutively. For example, Figure 5-1 shows the arrangement of data points for an image line of 1-byte data containing three spectral bands.
Figure 5-1. An image line may contain multispectral data stored in bands arranged consecutively.
|
For more information about manipulating data at the byte level, see the section titled Conversion utilities in Chapter 4, McIDAS-X Utilities. |
Below is a code fragment showing the steps for reading an image object data block. Note that the mcaget call must precede a call to mcalin .
If processing the image requires the line prefix, the application must get the prefix using the mcapfx function. The line prefix is the set of information that may precede the data on an image line; its maximum size is 1000 bytes. As shown in Figure 5-2, an image line consists of an optional line prefix and the actual data values.
Figure 5-2. Each image line contains an optional line prefix and data values.
The line prefix is divided into four regions:
The line prefix is the same length for each line in an image. The sensor source determines the size and content of the line prefix and which of its four regions are present.
As shown in the sample code below, the call to mcapfx must occur immediately after the call to mcalin so that the prefix and data are for the same image line.
The mcanav function reads an image object's navigation block. At any point after the connection is opened by mcaget , the application may retrieve the navigation block using the handle returned by the preceding mcaget call. See the code segment below. To eliminate the risk of buffer overflow, dimension the nav_buffer to 4096 bytes.
The mcacal function reads the calibration block of an image object. The call to mcacal can occur any time after the connection is opened by the mcaget call. The handle returned by mcaget is passed to mcacal , which returns the associated calibration block. To eliminate the risk of buffer overflow, dimension the cal_buffer to 40,000 bytes.
To read the comment block, use mcacrd and the handle returned by mcaget . The mcacrd function returns the entire comment block to the application. The call to mcacrd can occur only after the calls to mcalin are done.
Below is an example using mcaget , mcalin and mcacrd to read a comment block. To eliminate the risk of buffer overflow, dimension the comment_buffer to 40,000 bytes.
To write an image object to an image dataset, the application must identify a dataset and position, and open a connection with the server. Use the API functions below to write image objects to a dataset.
|
Function |
Description |
|---|---|
|
opens a connection to write the directory, navigation and calibration blocks of an image object |
|
The request to open a connection is performed by the mcaput function, which requires the following:
Because mcaput does not return an object handle, only one image object is written at a time.
The only valid selection clause for writing image objects is POS, which defines the location of the image object in the dataset. You must specify this clause or the request to open a connection will fail.
Once the connection is open, the server expects to transfer the number of bytes defined by the entries in the directory block. Transferring too few or too many bytes results in an error. All data block write transactions are performed by mcaout , which is called as many times as necessary to transfer the bytes. The mcaout function has only one argument, which is an array of data points to write to the data block. The call to mcaput must occur prior to mcaout .
The comment block is written after the last byte of the data block is transferred. The number of comment entries is defined in word 64 of the directory block. If this entry is nonzero, the specified number of 80-byte entries is transferred. The mcacou function transfers the comment block. It has only one argument, which is an array holding the entire comment block.
The sample code fragment below writes image objects to a dataset. For another mcaput code example, see the function imgcopy.f
c --- initialize the directory block call zeros(directory_block, 64) c --- create a comment card call getday( day ) call gettim( time ) cday = cfu( day ) ctime = cfu( time ) comment = cday(1:5)//' '//ctime(1:6)//' This is a comment ' ncard = ( len_trim(comment) / 80 ) + 1 c --- fill the essential directory block entries directory_block(2) = 4 ! version directory_block(3) = sss ! satellite number directory_block(4) = jday ! Julian day of image directory_block(5) = time ! nominal start time of image directory_block(6) = start_line ! starting image line number directory_block(7) = start_elem ! starting image element number directory_block(9) = num_lines ! number of lines of image data directory_block(10)= num_elems ! number of data points/line directory_block(11)= num_bytes ! number of bytes/data element directory_block(12)= line_res ! line resolution directory_block(13)= elem_res ! element resolution directory_block(14)= num_bands ! number of bands directory_block(19)= 2**(band-1) ! band map call movcw(memo, directory_block(25)) ! memo field directory_block(34)= data_offset ! byte offset to the data block directory_block(35)= nav_offset ! byte offset to the nav block directory_block(49)= doc_length ! length of prefix doc section directory_block(50)= cal_length ! length of prefix cal section directory_block(51)= lev_length ! length of prefix lev section directory_block(52)= lit( stype ) ! sensor source type directory_block(53)= lit( ctype ) ! calibration type directory_block(63)= cal_offset ! byte offset to the cal block directory_block(64)= ncard ! number of comment cards c --- initialize the navigation block call zeros(nav_block, nav_size) c --- fill the navigation block entries c Note: "navigation_params" is an array of navigation parameters c that describes the geo-location of the elements of the c image object. do 10 i = 1,nav_size nav_block(i) = navigation_params(i) 10 continue c --- initialize the calibration block call zeros(cal_block, cal_size) c --- fill the calibration block entries c Note: "calibration_parms" is an array of calibration parameters c that transforms the data elements to physical units. do 20 i = 1,cal_size cal_block(i) = calibration_params(i) 20 continue c --- fill the selection array nselect = 1 selects(nselect) = 'POS '//cfu(position) c --- open a connection to write the image object if( mcaput( image, nselect, selects, directory_block, nav_block, & cal_block).ne.0 ) then call edest('Unable to initialize image ='//image,0) return endif c --- loop to write image lines to the image object do 100 line = 1,num_lines c --- pack the data array c "data_array" is a (num_lines) by (num_elems) array of data c elements each of which is (num_bytes) long. The elements c represent data for (band) from the sensor numbered (sss) c on (jday) at (time). "data_buffer" is a one-dimension array c sized to (num_elems) c Note: this assumes a 4 byte to 1 byte compression of the data. call pack( num_elems, data_array(line, 1), data_buffer) c --- write a line of data to the image object if( mcaout( data_buffer ).ne.0 ) then call edest('failed to write image line=',line) return endif 100 continue c --- write the comment block if( mcacou( comment ).ne.0 ) then call edest('failed to write comment block',0) return endif
You can delete image objects from a dataset using mcadel, as long as you have write permission. The mcadel function has one valid selection clause, SUBSET, which you must specify during the connection phase of the transaction or the server request will fail. The code fragment below deletes the image objects at locations pos1 through pos2 from the dataset.
McIDAS-X grids are typically composed of atmospheric and oceanographic data, which are produced by numerical models or derived from observational data using an objective analysis scheme. Grids may contain any data that can be represented in a two-dimensional matrix.
Because grid data and image data can both be represented in two-dimensional matrices, it's important to know how they're different. The attributes that distinguish them are listed in the table below.
|
Attribute |
Grid data |
Image data |
|---|---|---|
Grid data has five unique attributes: valid time, level, parameter, origin and navigation. Each is described below.
Grid data is often used to store output from numerical models that simulate the atmosphere and ocean. Because these models predict what the atmosphere or ocean will look like at some time in the future, grids must contain two different time attributes when representing model forecasts:
For example, if a model run on 17 January at 00 UTC generates forecast fields valid 36 hours later, the primary time attribute is 17 January at 00 UTC and the secondary time attribute is 18 January at 12 UTC. If the grid does not represent a model forecast, no secondary time attribute is needed.
Grids store information that represents data at some vertical location in the atmosphere. The level can be the constant height above the surface, such as 100 meters or 32000 feet; the constant pressure surface, such as 500 millibars; or some other meteorological surface, such as the level of the tropopause or isentropic surface.
The grid parameter is the data type the grid represents. Parameters can be any field that can be stored as a number, such as temperature, wind speed or dew point. Grids must also contain the stored units used in the grid.
A grid's origin refers to the process that created the grid. For example, if a grid is created by a forecast model, the origin is the name of the model, such as ECMWF or GFS. If the grid is created from another objective analysis process, that name can appear as the origin.
Navigation is the process of determining the latitude and longitude location of data points on a grid. This information is needed when colocating data points from a grid with another data type. For example, you can use this information to contour the gridded field on a satellite image displayed on the McIDAS-X Image Window.
McIDAS-X currently recognizes these grid projection formats:
The diagrams below show how a three-dimensional Earth is represented on a two-dimensional surface for the pseudo-Mercator, polar stereographic and Lambert conformal projection formats, and the orientation of the data points on these projections. A sixth grid projection format is also available, but since it has no navigation, the data points can't be converted to planetary coordinates.
|
For more information about grid projections and McIDAS-X navigation, see the section titled McIDAS navigation later in this chapter. |
Figure 5-3. Pseudo-Mercator projection.
Figure 5-4. Polar stereographic projection .
Figure 5-5. Lambert conformal secant cone (left) and tangent cone (right) projections
Gridded data is two-dimensional data representing an atmospheric or oceanic parameter along an evenly spaced matrix. For the matrix to be useful, ancillary information about the grid must also be known. This ancillary information, along with the gridded data, is collectively called a grid object . Grid objects in McIDAS-X contain two blocks of information.
The API functions and procedures for reading and writing grid objects are described below.
Most applications for reading grid objects will do one of the following:
The table below lists alphabetically the McIDAS-X library functions for performing these tasks.
These functions are described below along with sample code, following an explanation of the selection conditions for requesting grid objects.
Applications use selection clauses to restrict the information sent from the server to the client. For example, selection clauses can restrict the search based on grid time attributes, the parameter type and level, or the process that generated the grid. Below is a list of the valid grid selection clauses as sent to the server. Additional information for each follows.
|
Selection clause |
Description |
|---|---|
DAY --use this clause to identify a list of primary days that the grids represent. For model forecast grids, these values will be the day the model was initialized. Otherwise, DAY is the Julian day the data represents.
DRANGE --use this clause to identify a range of primary days that the grid represents. Enter the beginning and ending day numbers and the increment between days in the range. The increment default is one day.
FDAY --use this clause to identify the secondary day attribute for requested forecast grids. For example, if you want only the forecast fields valid for day 1997017, specify FDAY=1997017.
FHOUR --use this clause to identify a list of secondary forecast hours that the grids represent. For example if you want only the 12-, 24- and 48-hour forecasts from a model run, specify FHOUR=12 24 48.
FRANGE --use this clause to identify a range of forecast hours that the grids represent. Enter the beginning and ending forecast hours and the increment between hours in the range. The increment default is one hour.
FTIME --use this clause to identify the secondary time attribute for the grids requested. Use this field with the FDAY clause to isolate grids that are valid at a particular time. The format for the arguments is hhmmss . For example, to request all grids valid at 12 UTC on day 96017, specify FDAY=96017 FTIME=120000.
GRIB --use this clause to select grids with the GRIB values specified for the geographic region, the parameter, the model, or the level. This is a useful way to find grids when other sort conditions provide insufficient information to differentiate between grids.
GRID --use this clause to access grids based on their position in specific grid files. Since this is an artifact of previous McIDAS-X API functions, avoid using this clause if a practical alternative exists.
LEV --use this clause to identify a list (not a range) of levels in the atmosphere or ocean that this data represents. This field is typically filled with height in millibars or words such as SFC, MSL or TRO. To retrieve data for several levels, enumerate them individually.
NUM --use this clause to specify the maximum number of grids returned from the server. The default is one grid. To receive all grids matching the selection conditions, use NUM=ALL.
OUT --use this clause to receive only the grid header. To get the entire grid header, specify ALL (default). To get only the grid file header, use FILE.
PARM --use this clause to specify a list of parameter types to retrieve from the server. To retrieve temperature and height fields, enter PARM=T Z.
POS --use this clause to specify a grid file in a dataset. This is a relative position based on the dataset description. For example, to request grid file 5010 from a dataset that contains grid files 5001 to 5100, specify POS 10.
PRO --use this clause to specify a projection type for a grid. The valid entries are LAMB, CONF and MERC.
SRC --use this clause to specify a list of grid source names to retrieve from the server. This is usually the name of the model or process that generated the grid, such as ETA, NGM or MDX.
TIME --use this clause to identify a list of primary times that the grids represent. For model forecast grids, these values will be the time of day that the model was initialized. The format for the arguments is hhmmss .
TRANGE --use this clause to identify a range of primary times that the grid represents. Enter the beginning and ending times and the increment between times in the range. The increment default is one hour. The format for the arguments is hhmmss .
You can use the m0gsort function with any application-level program to retrieve command line keyword parameters and translate them into equivalent selection clauses. The table below lists the keywords for accessing grid objects and their equivalent selection clauses. You can also set a flag in m0gsort to disable a request to contain multiple grid selection matches.
|
Command line keyword |
Equivalent selection clause |
m0gsort restrictions |
|---|---|---|
|
use LAST to get the last grid in the dataset and position; when specified, all other selection conditions are ignored |
||
In ADDE, a client may request only grid headers from a server. This allows the client to sample information on a server without transferring large amounts of grid data that may not be necessary for an application. For example, the McIDAS-X GRDLIST command reads only grid headers. The ADDE interface to the grid directory is through mcgdir , which opens a connection based on a set of selection clauses for a given dataset name.
Once mcgdir opens the connection, the application makes repeated calls to mcgfdrd and mcgdrd until all grid and grid file headers are retrieved. The mcgfdrd function is called first to retrieve the grid file header, then mcgdrd is called until it can't find any more grids in the grid file. Then mcgfdrd is called again and the loop continues, as shown below.
character*32 selects(5) character*32 dataset integer grid_header(64) integer file_header(64) integer nselects integer error_flag integer status c--- assign the dataset name dataset = 'RTGRIDS/ALL' c--- assign the selection conditions to retrieve six grids from c--- the dataset RTGRIDS/ALL that are from the ETA, NGM or c--- MRF model with a primary day of either 96017 or 96019 selects(1) = 'DAY=96017 96019' selects(2) = 'SRC=ETA NGM MRF' selects(3) = 'NUM=6' nselects = 3 c--- set an error flag to print a message if an error occurs error_flag = 1 c--- open the connection to the server status = mcgdir(dataset, nselect, selects, error_flag) if (status .lt. 0)then return endif c--- every time statement 100 is reached, try to read a new grid c--- file header 100 continue hdread = mcgfdrd(file_header) c--- if you have successfully read the grid file header if (hdread .eq. 0)then c--- every time statement 200 is reached, try to read a new grid header 200 continue gread = mcgdrd(grid_header) c--- if you have successfully read the grid header if (gread .eq. 0)then {process grid header here} c--- see if there are any more grid headers to read goto 200 c--- if you have read the last grid from this file, go see if there c--- are any more grid files to read from elseif (gread .eq. 1)then goto 100 c--- if there was a problem reading the grid header elseif (gread .lt. 0)then return endif elseif (hdread .lt. 0)then call sdest(`Unable to read grid file header',0) endif c--- if you make it to here, hread has returned the value 1, c--- which means the server has finished sending data
The grid header contains a list of ancillary information about the grid. The entries in the grid header are described in the table below.
|
Header Word |
Description |
|---|---|
|
total size; rows * columns (not to exceed the value of MAXGRIDPT in gridparm.inc) |
|
|
value of the vertical level
|
|
|
gridded variable type:
|
|
|
used if the grid parameter is a time difference or time average, hhmmss |
|
|
used if the grid parameter is a level difference or level average; values are the same as Word 9 |
|
|
grid origin; identifies the type of program that generated the grid data |
|
|
grid projection type:
|
|
|
varies, depending on the grid type; see the GRIDnnnn data file in Chapter 6 for more information |
|
|
reserved; filled only if the grid was created by the McIDAS-XCD GRIB decoder |
|
In ADDE, a client can request an entire grid object, including the grid header and data block from the server. For example, the McIDAS-X GRDDISP and GRDCOPY commands request entire grid objects.
The mcgget function opens a connection to read the data block of a grid object. It passes an application's selection conditions for requesting grid objects from the client to the server. The return status from mcgget indicates if the application's request can be fulfilled.
The mcgget function also allows the application to separately specify the units and format of the data returned. These are not part of the selection conditions because they are required. Units may be any unit identifier valid for the data type being retrieved. The format parameter can be either I4 for integer or R4 for real number.
If a grid in the dataset satisfies the client request, a connection is established between the server and the client and the transaction proceeds. The requested grid objects are read with either the mcgridf or mcgridc function; mcgridf reads the column-major (Fortran) format, while mcgridc reads the row-major (C) format.
Unlike image objects, which require multiple calls to mcalin to get an entire object, mcgridf and mcgridc return an entire grid object with each call. Below is a sample code fragment demonstrating a mcgget/mcgridf calling pair. Note that the mcgget call must occur before mcgridf .
Writing a grid object to a grid dataset has two restrictions:
Use the API functions below to write grid objects to a dataset.
|
Function |
Description |
|---|---|
|
writes a grid object stored in column-major format to a file |
|
The request to open a connection for writing a grid object is performed by the function mcgput , which requires the following:
When writing a grid object, the application may create and initialize a file in the destination dataset using the selection conditions below.
Once the connection is open, the server expects to transfer the number of bytes specified by the entries in the grid header. Transferring too few or too many bytes will result in a error. The mcgoutf and mcgoutc functions send the grid objects to the server. mcgoutf assumes the data block is stored in column-major order; mcgoutc assumes the data block is stored in row-major order. You must call mcgput before calling mcgoutf or mcgoutc , as shown in the sample code fragment below.
integer grid1(MAXGRIDPT), grid2(MAXGRIDPT) integer grid_h1(64) grid_h2(64) character*48 selects(10) character*48 dataset c--- write to the dataset LOCAL/GRIDS dataset = 'LOCAL/GRIDS' c--- set up the selection conditions. We will write 2 grid object c--- to dataset position 3 on the server. selects(1) = 'POS=3' selects(2) = 'NUM=2' selects(3) = 'MAX=100' nselects = 3 c--- initialize the grid headers. Not all fields are shown c--- in this example. The first grid object will contain 300 c--- mb height fields, the second will contain 250 mb height c--- fields c--- assign the size of the first grid object grid_h1(1) = 500 grid_h1(2) = 20 grid_h1(3) = 25 c--- assign 300 mb level grid_h1(10) = 300 : : c--- copy the first grid header to the second call movw(64, grid_h1, grid_h2) c--- change the level of the second grid header to 250 mb grid_h2(10) = 250 c--- assign the total number of bytes that will be transmitted. c--- this number will be the total number of data points in c--- each of the data blocks, plus the size of the grid header c--- of each of the objects. We multiply the total by 4 because c--- the storage format of the grid object is whole words and c--- there are 4 bytes per word. ngrids=2 total_bytes = 8+(grid_h1(1)+grid_h2(1)+(ngrids*64)+ngrids)*4 c--- send the request to write data to the server ok = mcgput(dataset, nselects, selects, 1, total_bytes) if (ok .lt. 0)then return endif c--- write the first grid object ok = mcgoutf(grid1, grid_h1) if (ok .lt. 0)then return endif c--- write the second grid object ok = mcgoutf(grid2, grid_h2) if (ok .lt. 0)then return endif
In the sample code above, you will notice that the last parameter in the mcgput function is the total number of bytes to transfer from the client to the server. You must provide this number using the equation below.
total bytes = 8 + 4 * ((number of grids * 64) + (number of data points ) + (number of grids ))
For example, to transfer two grids to a server where the first grid has 200 data points and the second grid has 300 data points, the total number of bytes to transfer is 2528, as shown below.
McIDAS-X point data is typically composed of atmospheric and oceanographic data occurring at irregularly spaced locations on the Earth or vertically within the atmosphere or ocean. This type of data storage is most often used with station observations such as synoptic, RAOB or ship reports.
Point data has two unique attributes:
Each of these attributes is described below along with the limitations imposed on point data.
Conceptually, you can think of point data as a spreadsheet with each cell containing a predefined number of data values. Each cell contains data for a specific location at a given instant in time. For example, one cell may contain all the mandatory level RAOB data from Green Bay, Wisconsin, at 12 UTC on 17 January 1996.
Unlike grids, point objects may contain a combination of data types and units for different elements within a cell. For example, in a point cell of surface hourly information, you can store the following:
When accessing point data, the MD (Meteorological Data) file structure has four limitations, which are explained below. These limitations pertain only to the MD file structure and are not limitations of the ADDE point object subsystem. Point servers for other formats, netCDF for example, use configuration files to map the parameter names in the file to McIDAS-compatible names.
Each cell is limited to 400 elements. For most observational data, this isn't problematic. The exception occurs with observations containing data at several levels of the atmosphere. For example, you can't store all the information for an upper air observation reporting values for eight different parameters at 50 levels of the atmosphere in one cell. Although the cell can accommodate the 400 values (8 x 50), it won't have enough space for the time and geographic location of the observation, which are also provided.
Element names and units are limited to four characters, which can be restricting when designating parameter units, especially derived parameters.
Character string elements are limited to four characters, which can be limiting for any type of alphanumeric parameter, such as station ID or country. You can bypass this restriction by using several parameters strung together to represent strings. McIDAS-XCD uses this method when representing five-character IDs associated with ship reports.
McIDAS-X point data typically represents data occurring at irregularly spaced locations on the Earth. For this data to be useful, ancillary information about the data must also be known. This ancillary information, combined with the actual point data values, is collectively called the
point object .
Each point object in McIDAS-X contains five blocks of information.
The API functions and the procedures for reading point objects are described below.
Most applications for reading point objects will request one of the following:
The McIDAS-X command PTLIST is an example of an application that reads point objects. The table below lists alphabetically the McIDAS-X library functions for acquiring point data.
|
Function |
Description |
|---|---|
|
gets the selection parameters from the command line and adds them to a selection array used by m0ptget |
These functions are described below along with an explanation of the selection conditions for requesting point objects.
Applications use selection clauses to restrict the information sent from the server to the client. You can tell the server to return only fields that fall within certain thresholds. These selection limitations may include a list of stations, a time range, or a level in the atmosphere. Below is a list of the valid point selection clauses. Additional information for each follows.
|
Selection clause |
Description |
|---|---|
|
position number in a dataset to retrieve data from (default=ALL if day is specified; default=0 if day is not specified) |
|
MAX-- use this clause to specify the maximum number of point objects returned from the server. To receive all point data matching the selection conditions, use MAX=ALL.
POS-- use this clause to identify a specific point file in a dataset. This is a relative position based on the dataset description. For example, to request point file 5010 from a dataset containing point files 5001 to 5100, specify POS=10. To read from all the files in a dataset, specify POS=ALL.
SELECT-- use this clause to identify the selection conditions for limiting the objects returned to the client. The syntax of this clause varies, depending on the request. To include multiple conditions, enclose each clause in single quotes. For example, to limit the list of stations in a surface hourly observation to include only Madison and Milwaukee, you can use: 'ID KMSN, KMKE'. To limit a selection to include only those stations with a temperature between 0° and 32° Fahrenheit, include: 'T 0 TO 32 F'.
Since piecing the selection conditions together can be a difficult task, the m0psort function will build the appropriate SELECT string for you. Use m0psort with any application-level program to retrieve command line keyword parameters and translate them into equivalent selection clauses. For example, if you enter a point command with the following arguments:
m0psort will return the following information, which can then be passed to the server via m0ptget.
Once the selection conditions are made, the m0ptget function opens a connection to read point data from the server. The calling sequence for m0ptget allows the client to access the data in one of two modes.
In the first mode, the server returns all the parameters for a given data type matching the selection conditions, such as all the decoded information from a METAR report for a given station. If you use the PTLIST command without specifying the PARAM keyword in this mode, the client may not know the valid parameters and units for the data type until a successful return from m0ptget .
In the second mode, the client knows the list of parameters to request and the units they can be returned in, such as temperature in Celsius and wind speed in knots.
The API for m0ptget contains the field asknparm for input and output. If this field is zero, the first mode of data acquisition is assumed and the client retrieves all data associated with this data type. If a list of parameters is specified from the client, the second mode of data acquisition occurs and asknparm will contain the number of parameters to return.
When accessing point data in the second mode, you must supply m0ptget with a list of parameters and units to retrieve. If you adhere to the McIDAS-X command line syntax PARAM=parameters [units ], you can use the m0ptparm function to build the parameter and unit list.
Upon successful return from m0ptget , all elements of the point object are returned by the server, except the actual data block. If the client requests all parameters, the parameter and unit block are returned. The data form is also returned, indicating the type of data each parameter is stored in. For character strings, the form is C#; for integer values, it is I#; and for floating point values, it is F# , where # is the number of bytes for this parameter. The scaling factor is also returned for floating point numbers.
If the request can be fulfilled by the server, the parameter, unit, scale and form blocks are filled accordingly, regardless of which mode the client uses for m0ptget . The example in the next section demonstrates the second mode of data retrieval.
If the point request for a given dataset can be fulfilled, a connection is established between the server and the client and the transaction proceeds. The m0ptrd function reads the point data. It is called continuously until all the data is read. Because each call to m0ptrd may yield character strings, integers and floating point numbers, the McIDAS-X library contains a group of functions for extracting these mixed data type values. You can call the m0ptbufinit function after each m0ptget call to initialize the application environment to more easily extract data from the buffer.
The example below demonstrates point data acquisition, using the m0ptbufc , m0ptbuff and m0ptbufi functions to extract characters, floating point, and integer values from the buffer filled by m0ptrd .
subroutine main0 implicit none include 'ptparm.inc' integer MAXBYTE ! max bytes the buffer may contain parameter (MAXBYTE = MAXNUMPARM * 4) character*24 dataset ! adde dataset name character*80 select(2) ! list of selection conditions integer nselect ! number of selection conditions integer nparams ! number of parameters to return character*1 quote ! single quote character*4 params(MAXNUMPARM) ! list of parameters character*4 units(MAXNUMPARM) ! list of units character*4 form(MAXNUMPARM) ! list of return forms integer scales(MAXNUMPARM) ! list of scales integer buffer(MAXNUMPARM) ! data buffer integer ok ! function return value double precision temperature ! temperature extracted from the ! buffer double precision cloudheight ! cloud height extracted from the ! buffer character*4 id ! station id extracted from the ! buffer character*12 ct ! temperature character*12 czcl ! cloud height character*12 chms ! ob time integer obtime ! observation time extracted ! from the buffer c--- external library routines integer m0ptbufinit integer m0ptbufc integer m0ptbufi integer m0ptbuff integer m0ptget integer m0ptrd data units/MAXNUMPARM * ' '/ quote = char (39) c--- assign the dataset dataset = 'RTPTSRC/SFCHOURLY' c--- set up the selection conditions to retrieve the station id, c--- the temperature in Celsius, the observation time and the c--- first non-ceiling cloud height for five SFCHOURLY stations c--- in Wisconsin with a temperature between 60 and 70 Fahrenheit c--- for day 1996274 params(1) = 'ID' params(2) = 'T' units(2) = 'C' params(3) = 'HMS' params(4) = 'ZCL1' nparams = 4 select(1) = 'MAX=5' c--- the assignment below is usually more easily performed by the c--- function m0psort. it is demonstrated in this example in a fully c--- expanded manner to show the most basic access method. select(2) = 'SELECT=' & //quote//'ST WI'//quote//' ' & //quote//'T 60 TO 70 F'//quote//' ' & //quote//'DAY 1999274'//quote c--- select(2) will contain the following. note that the single c--- quote characters are an important element of the string: c--- SELECT='ST WI' 'T 60 TO 70 F' 'DAY 1999274' nselect = 2 c--- make the request to the server ok = m0ptget (dataset, nselect, select, nparams, params, & units, form, scales, MAXBYTE, 1) if (ok .lt. 0)then goto 999 endif c--- upon a successful return from m0ptget, the arrays will contain c--- the following: c--- location params units form scales c--- 1 ID CHAR C4 0 c--- 2 T C F4 2 c--- 3 HMS HMS I4 0 c--- 4 ZCL1 FT F4 -2 c--- initialize the system so m0ptbufc/f/i can be used to c--- extract the buffer values ok = m0ptbufinit (nparams, form, scales) if (ok .lt. 0)then goto 999 endif c--- continuously call m0ptrd until it indicates there are no more c--- point data blocks to return. call sdest('ID T[C] Time ZCL1[FT]',0) 100 continue ok = m0ptrd (buffer) if (ok .lt. 0)then goto 999 endif c--- if data was found, process the data if (ok .eq. 0)then ct = ' ' chms = ' ' czcl = ' ' c--- extract the station id ok = m0ptbufc (1, buffer, id) c--- extract the temperature ok = m0ptbuff (2, buffer, temperature) if (ok .eq. 0)then write(ct,FMT='(f9.1)')temperature endif c--- extract the observation time ok = m0ptbufi (3, buffer, obtime) if (ok .eq. 0)then write(chms,FMT='(i6)')obtime endif c--- extract the cloud height ok = m0ptbuff (4, buffer, cloudheight) if (ok .eq. 0)then write(czcl,FMT='(f7.1)')cloudheight endif write(cline, & FMT='(a4,2x,a5,1x,a6,1x,a6)') & id,ct(6:),chms,czcl call sdest(cline,0) goto 100 endif 999 continue call edest('done',0) return end
Occasionally, an application may need to access the file header associated with point data without accessing the data itself. The McIDAS-X ADDE command PTLIST with the FORM=FILE option is an example of such an application.
The m0pthdr function opens a connection to read point data file headers from the server based on the selection conditions shown in the table below.
|
Selection clause |
Description |
|---|---|
Upon a successful return from m0pthdr , the m0ptrdhdr function is repeatedly called until all the data is returned. Below is a sample code fragment demonstrating the use of m0pthdr and m0ptrdhdr .
The contents of a header returned from m0ptrdhdr are described in the table below.
|
Word |
Description |
|---|---|
|
user record, MD coordinates (0,0); not described by the schema; use for storing arbitrary information |
|
Words 39-42 are filled in only during the production of real-time files. When real-time file data is copied to other MD files, words 39-42 in the destination files are set to null.
Text data stores a variety of information for the McIDAS-X user, such as:
McIDAS-X has two types of text data: flat-file text and general weather text. Both are described in this section along with:
Flat-file text is the simplest text data available in McIDAS-X. It is used most often on the server machine to convey administrative or configuration information to the user. The McIDAS-X READ command accesses this text.
The data is stored on the server machine as a simple ASCII file that can be manipulated with any file editor. It is delivered to the client one line at a time. There is no practical limit to the length of an individual line.
The table below lists the McIDAS-X library functions for flat file text.
|
Function |
Description |
|---|---|
|
opens a connection to read a flat file and prints the file contents to the McIDAS-X Text Window |
Because there is a one-to-one relationship between the ADDE dataset name created on the server and the file being served, no selection conditions are needed to access most flat files. Only the dataset name is required.
The ADDE interface to flat-file text is through the function M0textgt . It has one selection condition, FILE=, which can be used to access files that don't have specific dataset names in ADDE. The WXTLIST command with the DIR option is the only command that uses this selection condition.
Once the connection is opened with M0textgt , the M0txtread function is called continuously until no more lines of text data are found. The code example below demonstrates the M0textgt/M0txtread function pair.
General weather text is composed of information compiled by weather agencies and distributed to the user community. It typically contains forecasts, public announcements, advisories and warnings. In McIDAS-X, this type of data is ingested and made available from a server running McIDAS-XCD and can be accessed by the user with the McIDAS-X command WXTLIST.
Because weather agencies provide this data, certain WMO format standards must be used during transmission. The sample text below demonstrates the WMO formatting standards used for text data transmission.
The first line of text is the WMO header line. In the example above, the first two characters, FP, contain the WMO product header. Most WMO product headers begin with one of the following characters:
|
WMO product header |
Description |
|---|---|
The second character of the WMO product header will vary. The third and fourth characters are usually the country code (US) of the bulletin. The country code is followed by a product number (5), which further specifies the type of bulletin being transmitted. The product number is followed by the ID of the station initiating the bulletin (KMKE). The final six digits are the day of the month and time of the bulletin (251641).
The second line of the bulletin is an AWIPS header inserted by the U.S. National Weather Service; it is optional. The first three characters are the product code (SFP). The next two or three characters contain the state/province or the station the report is valid for (WI).
The table below lists the McIDAS-X library functions for reading general weather text.
|
Function |
Description |
|---|---|
In ADDE, the M0wtxget function opens a connection to the server. This function takes a variety of selection conditions, allowing the user to limit the amount of data searched to fulfill a user request. The more information specified in the selection conditions, the faster the search will be. The valid selection conditions are described in the table below.
|
Selection clause |
Description |
Remarks |
|---|---|---|
|
if more than one match string is requested, all match strings requested must be found in each text block for a successful return |
||
|
minimum of two characters; wildcard characters are allowed; see the help for the WXTLIST command; can't be used with APRO= |
||
Once M0wtxget opens the connection, the application makes repeated calls to M0wtxread until all matching text products are received. Each call to M0wtxread returns the actual text along with a 64-byte header containing information about the text data. The components of the header are described in the table below. Note that all character strings are sent blank padded.
|
Bytes |
Description |
|---|---|
The sample code below demonstrates an ADDE general weather text request.
All geophysical data must have a location to be meaningful. McIDAS-X stores geophysical data in the form closest to its native spatial structure; for example, remotely sensed data is stored in the image coordinates of the scanner line and element, while NWP model output is stored as model grid points. Since earth coordinates (for example, latitude, longitude and height above mean sea level) are easiest to interpret, McIDAS-X uses a process called navigation to convert a dataset's native coordinate system to earth coordinates and back again.
Navigating datasets implies the ability to transform any dataset into the spatial coordinates of any other dataset. McIDAS-X' ability to integrate and display a variety of datasets is made possible by its navigation subsystems.
|
For more information about McIDAS-X coordinate systems, see Chapter 2, Learning the Basics . For additional information about map projections, see Map Projections -- A Working Manual , USGS Professional Paper 1395, U.S. Government Printing Office, Washington, 1987. |
To select and use the best McIDAS-X navigation subsystem for a given task, you need to understand the nature of navigation transforms. Navigation systems and transforms are described below, along with some navigation-specific terminology used throughout this section.
A navigation transform is a set of equations for converting a dataset's image or native coordinates to and from earth coordinates. Each dataset has its own navigation transform. However, navigation transforms can be grouped into classes or types . Within a type, the differences between transforms are quantitative only. Each type has a set of variables, or parameters, whose values uniquely specify an instance of that type. For example, the parameters for a tangent cone projection include the standard latitude and longitude, the scale factor, and the location of the pole relative to the projection surface's origin.
All navigated McIDAS-X data has an associated type and a full set of parameters defining an instance.
Grid navigation converts earth coordinates to and from the row and column locations in the grid structure. This subsystem is separate from image and frame navigation, and supports fewer projections These projection types are implemented in a single set of functions. To add other types, you would have to modify the functions and recompile all grid applications. Since only a single instance can be in use at one time, you must repeatedly reinitialize the subsystem to use multiple projections within an application.
For example, to compare two grids, A and B, you initialize grid navigation for A, compute the earth location of a grid point a in A, reinitialize the subsystem using B's parameters, and transform the earth location of a to a grid in B. This process is then repeated for every grid point being compared. Most applications that process or use grids use grid navigation.
|
For more information about grid projection formats, see the Grid data section earlier in this chapter. |
Image navigation converts image coordinates to and from earth coordinates. Each image navigation type is implemented in its own source file or module , and each navigable image dataset has an associated navigation block that identifies the type and contains a complete set of parameters needed to define an instance of that type.
Each module implements the same function names and calling sequence. Since McIDAS-X uses the type identifier to select the correct module, you can add new types without modifying existing modules. In addition, up to three instances of navigation, of the same or different types, can be active simultaneously in a single application through different slots. This requires a one- or two-step initialization process before use and is typically used by applications that work directly with image datasets.
Frame navigation converts frame coordinates to and from earth coordinates. Frame navigation uses the image navigation subsystem but through a simpler API. It is sufficient for most applications that just add information to an existing display.
This section contains some navigation-specific terms that may be unfamiliar to you. They are defined below.
The McIDAS-X library provides Application Program Interfaces for geographic calculations as well as grid, image and frame navigation. These APIs are described below.
Many calculations just require a change from one form of earth or planetary coordinate to another; for example, changing from a terrestrial (earth-relative) to celestial (star-relative) coordinate system, or converting a position from a Cartesian to spherical representation. The table below lists the functions for performing these tasks.
|
Functions |
Description |
|---|---|
The table below shows the grid navigation API. You will use the grddef function with the header to initialize the subsystem, then call ijll and llij to transform grid coordinates to and from earth coordinates, respectively.
|
Function |
Description |
|---|---|
No API functions currently exist for creating grid headers for a projection. When creating McIDAS-X grid files, you must determine the projection to use and then insert the projection type and parameters into the grid header. The grid header format is described in the GRIDnnnn data structure in Chapter 6.
Use the two functions below to determine the earth coordinates of a pixel or the pixel location of a point in earth coordinates.
|
Function |
Description |
|---|---|
The only prerequisite for these functions is that the frame has navigation associated with it. The McIDAS-X IMGDISP command creates frame navigation when it displays an image dataset. The McIDAS-X MAP, PTDISP and GRDDISP commands, which can define a map base before generating a display, also establish frame navigation.
The two examples below assume that frame navigation information exists. The first code fragment uses the illtv function to convert an earth location (zlat,zlon ) to a frame location (itvlin,itvele ), as is done in the McIDAS-X PC command.
zlat = 43.13 zlon = 89.35 if(illtv(iframe,gnum,zlat,zlon,ilin,iele,itvlin,itvele).ne.0)then call edest('unable to navigate point',0) return endif |
The second code fragment uses itvll to determine the earth location of the cursor, as is done in the McIDAS-X E command.
if( mcmoubtn( 0, left, right, itvlin, itvele ) ) then C // handle cursor-positioning error here end if frame = mcgetimageframenumber() rc = itvll(frame,itvlin,itvele,ilin,iele,rlat,rlon,iscene) |
Use the appropriate image navigation API functions from the table below if your application must perform one of these tasks:
To use image navigation services, the application must do the following:
The example below, from remap.pgm, illustrates the initialization and use of two instances (slots) of navigation at the same time.
The lines below assume that navigation blocks were already obtained from the source (snav ) and destination (dnav ) areas. To link the appropriate navigation module to a slot and initialize it, use nvprep . The lines below initialize slot 1 for the source and slot 2 for the destination. The nv1ini and nv2ini calls set the earth coordinate form to Cartesian for both slots.
The lines below map image coordinate (dline,delement) in the destination image to image coordinate (sline,selement) in the source image by converting from destination image coordinates to earth coordinates, then from earth coordinates to source image coordinates.
When the navigation block source is a local McIDAS-X area format image file or a frame whose number is known and you only intend to use slot 1, you can use nvset. It reads the navigation block out of an area or frame and calls nvprep.
In addition to the commonly used conversions between earth and image coordinates available through nvxsae and nvxeas, other special services are available through the nvxopt function. Some of these services and the modules that support them are listed in the table below.
|
Special services |
Modules that support them |
|---|---|
|
nvxdmsp, nvxgoes, nvxgraf, nvxgvar, nvxlamb, nvxmerc, nvxmsat, nvxps, nvxradr, nvxrect, nvxsin, nvxtiro |
|
The sample code below shows the use of nvxopt. When using nvxopt, you need to be aware that the input/output can be different for the different satellite navigations, and you'll have to check the nvx sat.dlm file to check what is appropriate for input/output.
You can extend image navigation to include a new type and make it available to all McIDAS-X commands without changing any applications or existing navigation modules. This section describes the process for implementing a tangent-cone map projection. The complete source listing is in the next section titled Example navigation module code .
The navigation block must consist of an array of no more than 640 Fortran integers that are sufficient to uniquely specify a particular instance of the new navigation type. The first word of the block must be four ASCII characters that uniquely identify the type, stored in an integer variable. The block contents are obvious once the algorithm is fully specified and understood.
The navigation module name must be nvx type .dlm , where type is the same four characters specified in the first word of the block. The nvprep function uses the type from the block to link the application with the correct module.
The navigation module must implement the functions nvxini , nvxsae , nvxeas and nvxopt , with calling sequences following the documentation blocks in the tangent cone example. In particular, nvxini must support modes 1 (initialize) and 2 (set earth coordinate mode 'LL' (Latitude-longitude) or 'XYZ' (Cartesian)). The nvxopt function must provide, at minimum, an 'SPOS' (subpoint) service.
Between calls, the navigation module must remember the navigation parameters uniquely defining an instance. This is typically accomplished by storing these parameters, and quantities derived from them, in a named common block that is shared among all four API functions in the module.
Satellites and map projections cannot physically represent all earth locations. The limits vary from type to type and even from instance to instance. Applications, therefore, cannot be expected to validate inputs to navigation modules. The new module must be able to recognize inputs it cannot handle and return an error status, typically a negative integer.
A navigation algorithm is a set of equations for mapping McIDAS-X image coordinates to earth coordinates and vice versa. For satellite navigation types, this typically involves a lengthy expression in vector notation. Map projections are simpler. Some manipulation and derivation of additional relationships is often necessary to extend a published relationship into a full algorithm.
The example presented here will show the development of a northern hemisphere tangent cone projection algorithm. The geometry of the tangent cone is shown in Figure 5-6 below. The cone is tangent to the planet at one parallel of latitude and the planetary axis pierces the apex of the cone. A line between a point on the Earth's surface and the pole opposite the cone apex maps that point onto the cone. The imaginary cone is then cut along a meridian opposite an arbitrary standard longitude and flattened. Each point on the Earth's surface in colatitude ψ and longitude λ corresponds to a point (R,θ ) in polar coordinates on the flattened cone.
Figure 5-6. Tangent cone projection geometry .
A published form of the equations defining the location (R,θ) on the flattened cone in terms of an earth location (ψ,λ) is:
Also provided is an equation for the ratio of the distance between two closely spaced points on the cone to the distance between the corresponding points on the Earth's surface. As you can see from (1) and (2), the planetary radius, a the standard colatitude ψ0, and the standard longitude λ0 define an instance of the tangent cone. In equation 3, m is the map scale expressed as a unitless ratio of projection plane distance to map distance.This is not yet a complete algorithm ready for McIDAS-X implementation. The inverse relationship capable of converting earth location to projection coordinate (R,θ ) must be derived, the earth location (λ ,ψ ) related to its equivalent in McIDAS-X (LAT, LON), and the projection coordinate (R,θ ) related to McIDAS-X image coordinates. The inverse is:
and The relationships between McIDAS-X planetary coordinates and (λ ,ψ ) are: and Those between projection and image coordinates are: and Adding (6) - (9) to make a complete McIDAS-X navigation algorithm adds additional parameters to make the complete set ψ0, λ0, m, L0, and E0, where the latter specify the scale in km per pixel and the location of the pole of the projection surface in image coordinates.The last step before implementation is to examine the completed algorithm for limits and singularities; these will serve as the basis for input validation. The limits on the inputs are summarized below.
The contents of the navigation block follow directly from the parameters provided in the previous Algorithms section. The five parameters are:
A sixth parameter, planetary radius, could also be added but is not included in the examples. The parameters can be represented in different ways, such as standard colatitude in place of latitude. However, since the applications will create navigation blocks, you should use units and quantities convenient for McIDAS-X, as shown in the lines of code below. This code is taken from the documentation block at the beginning of nvxini and is shown in lines 0017 through 0022 in the complete source listing that follows.
All instances of 'TANC' navigation must be specified with these parameters.
Applications that create navigation blocks must provide values for a complete set of parameters. Sometimes they must be derived from user inputs using the navigation transformation itself. Below is a code fragment from an application (not the navigation module itself) that creates a tangent-cone navigation block. It assumes that the center latitude and longitude (clat and clon ), the standard latitude and longitude (slat and slon ), and the scale of standard latitude (sscale ), in km per pixel, are already provided. First, the earth locations are converted from McIDAS-X to projection form.
psi = PI/2.d0 - D2R*clat psi_0 = PI/2.d0 - D2R*slat lam = - D2R*clon lam_0 = - D2R*slon |
Next, the radius of the specified image coordinate center point (relative to the projection center or the pole) is computed.
radius = A * tan(psi_0) * $ (tan(psi/2.d0)/tan(psi_0/2.d0))**cos(psi_0) theta = cos(psi_0)*(lam-lam_0) |
Then the pole's image coordinate is computed using simple trigonometry.
call status = mcfsize(mcgetimageframenumber, nlins, neles) lin_c = dble(nlins) / 2.d0 ele_c = dble(neles) / 2.d0 lin_0 = lin_c - radius*cos(theta)/sscale ele_0 = ele_c - radius*sin(theta)/sscale |
Once the parameters are computed, the navigation block can be filled and written either to the image dataset:
or a McIDAS-X area format image file:
call araput(aranum,4*DIRSIZ,4*NVWDS,nvblk) |
If the algorithm and navigation block are complete and unambiguous, implementing a map projection as a McIDAS-X navigation module is generally straightforward. Your major design decision is what to store in the common block. At minimum, you should include the navigation parameters or quantities derived from them without loss of information. This is necessary for the module to remember the characteristics of the instance.
Because navigation modules are called many times, you should precompute additional quantities in the initialization routine nvxini and store them in common, especially if doing so eliminates trigonometric function calls in the other functions. Thus, the variables Coscl , Tancl , Tancl2 , and Mxtheta are included in the common block below. This code is located in lines 0099 through 0120 of the source listing that follows.
Also included in the block are physical constants and two flags; Init is set when nvxini runs successfully, and Latlon is set if the earth coordinate mode is latitude and longitude rather than Cartesian.
The overall structure of nvxini is a large if-block to handle the two options. Option 1 (lines 0144 - 0220 in the source listing) initializes the instance by range checking the navigation parameters and precomputing the quantities to retain in the common block. Option 2 (lines 0222 - 0233) sets the earth coordinate mode based on word param(1) .
The forward (image to Earth) navigation routine's one subtlety is input validation. Points within the split of the projection plane, which is shown in Figure 5-6, are not navigable.
Recognizing such points is possible only after the conversion from McIDAS-X image coordinates to projection coordinates is complete. The code below is from lines 0370 through 0384.
When range checking is done, the projection is inverted and the colatitude and longitude (in radians, east positive) is converted to latitude and longitude (in degrees, west positive); see lines 0389 - 0401. A call to nllxyz to generate Cartesian coordinate output is made if the Latlon flag was cleared by an earlier call to nvxini (line 0404). Inverse (earth-to-image) navigation in nvxeas follows a similar pattern, except the earth-coordinate option must be handled at the beginning (line 0544).
The nvxopt function performs special services. Since the contents of the argument vector depend on the option selected, the routine consists of a large if-block (lines 0716 - 0742) with one branch per recognized option and the necessary input validation done within each branch. As in nvxini , the code must be able to recognize options that it doesn't know and return an appropriate error status. All core McIDAS-X navigation modules have an SPOS option that returns the subpoint at a given time (for satellites) or the latitude and longitude of the center (maps). The example that follows implements only an SCAL to return the map scale factor at any latitude.
1. Saucier, W. J. 1983: Principles of meteorological analysis . Dover Publications, Inc., New York. 438 pp.
0001 C THIS IS SSEC PROPRIETARY SOFTWARE - ITS USE IS RESTRICTED. 0002 0003 C *** McIDAS Revision History *** 0004 C *** McIDAS Revision History *** 0005 0006 *$ Name: 0007 *$ nvxini - Initialize navigation for tangent cone projection 0008 *$ 0009 *$ Interface: 0010 *$ integer function 0011 *$ nvxint(integer option, integer param(*)) 0012 *$ 0013 *$ Input: 0014 *$ option - 1 to set or change projection parameters 0015 *$ option - 2 set output option 0016 *$ param - For option 1: 0017 *$ param( 1) = 'TANC' 0018 *$ param( 2) = image line of pole*10000 0019 *$ param( 3) = image element of pole*10000 0020 *$ param( 4) = km per pixel *10000 0021 *$ param( 6) = standard latitude *10000 0022 *$ param( 7) = standard longitude *10000 0023 *$ for option 2: 0024 *$ param( 1) = 'LL' or 'XYZ' 0025 *$ 0026 *$ Input and Output: 0027 *$ none 0028 *$ 0029 *$ Output: 0030 *$ none 0031 *$ 0032 *$ Return values: 0033 *$ 0 - success 0034 *$ -3 - invalid or inconsistent navigation parameters 0035 *$ -4 - invalid navigation parameter type 0036 *$ -5 - invalid nvxini() option 0037 *$ 0038 *$ Remarks: 0039 *$ Latitudes and longitudes are in degrees, West positive. 0040 *$ Projection parameters must be in the following ranges: 0041 *$ 0. < standard latitude < 90. 0042 *$ -180. <= standard longitude < 180. 0043 *$ 0. < scale 0044 *$ Accuracy may suffer near the standard latitude limits. 0045 *$ 0046 *$ The projection algorithm is adapted from that in 0047 *$ Saucier, W. J. 1989: Principles of meteorological analysis. 0048 *$ Dover Publications, Inc. 433 pp. 0049 *$ 0050 *$ Categories: 0051 *$ navigation 0052 0053 C // CODING CONVENTION note: function declarations and common 0054 C // block declarations are all capitalized to be recognizable 0055 C // to script 'convdlm;' this is necessary for a correct build 0056 C // in MCIDAS-X. For the same reason, avoid referring to 0057 C // function or common block names in uppercase elsewhere 0058 0059 INTEGER FUNCTION NVXINI(option,param) 0060 0061 0062 implicit NONE 0063 0064 0065 C // Interface variables (formal arguments) 0066 0067 integer option ! initialization option 0068 integer param(*) ! navigation parameters or 0069 C ! output coordinate type 0070 0071 C // Local variable definitions 0072 0073 character*4 navtyp ! codicil type 0074 character*4 outcoord ! output coordinate type 0075 real lat0 ! standard latitude 0076 character*80 cbuf ! text output buffer 0077 0078 0079 0080 C ///////////////////////////////////////////////////////// 0081 C Common block variables and declaration. 0082 0083 C ALL CODE BETWEEN THE '//////' SEPARATORS MUST BE 0084 C DUPLICATED EXACTLY IN EACH NAVIGATION ROUTINE 0085 0086 C (A more maintenance-safe version would use ENTRY points 0087 C rather than separate functions for the navigation API 0088 C but entry points cannot be processed by 'convdlm.') 0089 0090 C // Common block contents: projection parameters 0091 0092 real Lin0 ! image line of pole 0093 real Ele0 ! image element of pole 0094 real Scale ! km per unit image coordinate 0095 C ! (pixel) 0096 real Lon0 ! standard longitude 0097 real Colat0 ! standard colatitude 0098 0099 C // Common block contents: precomputed intermediate values 0100 0101 real Coscl ! cosine(Colat0) 0102 real Tancl ! tangent(Colat0) 0103 real Tancl2 ! tangent(Colat0/2) 0104 real Mxtheta ! limit of angle from std. lon 0105 C ! on projection surface 0106 0107 C // Common block contents: constants 0108 0109 real D2R ! degrees to radians factor 0110 real Pi 0111 real Badreal ! returned when navigation 0112 C ! cannot be done 0113 real Erad ! Earth radius 0114 logical Init ! initialized flag 0115 logical Latlon ! .TRUE. for lat/lon I/O 0116 0117 0118 COMMON/TANC/ Lin0, Ele0, Scale, Lon0, Colat0, 0119 & Coscl, Tancl, Tancl2, Mxtheta, 0120 & D2R, Pi, Badreal, Erad, Init, Latlon 0121 0122 C End of common block variables and declaration. 0123 C ///////////////////////////////////////////////////////// 0124 0125 0126 0127 C // Begin initialization process by setting values of constants. 0128 0129 Erad = 6370. ! This value of Erad is ok 0130 C ! for low precision nav 0131 C ! where a spherical Earth is 0132 C ! adequate (Saucier, p. 32) 0133 Pi = acos(-1.) 0134 D2R = Pi / 180. 0135 Badreal = -1.E10 ! obvious unreasonable value 0136 C ! for nav transform result 0137 0138 C // Process initialization options. Only one option, initialize 0139 C // navigation parameters, is supported in this demo version, 0140 C // but a 'hook' is left for an additional option to set the 0141 C // output coordinate to something other than lat/lon 0142 0143 if( option.eq.1 ) then 0144 0145 call DDEST('nvxini(tanc) option=1',0) 0146 0147 call movwc(param(1),navtyp) 0148 if( navtyp.eq.'TANC') then 0149 0150 C // Unpack tangent cone projection parameters 0151 0152 Lin0 = param(2) / 10000. 0153 Ele0 = param(3) / 10000. 0154 Scale = param(4) / 10000. 0155 lat0 = param(5) / 10000. 0156 Lon0 = param(6) / 10000. 0157 0158 write(cbuf,'('' nvxini: lat0, Lon0 '',2F12.4)') 0159 * lat0, Lon0 0160 call DDEST(cbuf,0) 0161 0162 C // apply range checking 0163 0164 if(Scale.le.0. ) then 0165 call DDEST('nvxini(tanc) scale is negative',0) 0166 Init = .FALSE. 0167 NVXINI = -3 0168 return 0169 end if 0170 0171 if(lat0.le.0. .or. lat0.ge.90. ) then 0172 call DDEST('nvxini(tanc) standard lat out of range',0) 0173 Init = .FALSE. 0174 NVXINI = -3 0175 return 0176 end if 0177 0178 if(Lon0.le.-180. .or. Lon0.gt.180. ) then 0179 call DDEST('nvxini(tanc) standard lon out of range',0) 0180 Init = .FALSE. 0181 NVXINI = -3 0182 return 0183 end if 0184 0185 C // convert degrees to radians and latitude to colatitude. 0186 C // Account for McIDAS longitude convention 0187 0188 Lon0 = -Lon0 * D2R 0189 Colat0 = Pi/2. - D2R*lat0 0190 0191 write(cbuf,'('' nvxini: Colat0, Lon_0 '',2F12.4)') 0192 * Colat0, Lon0 0193 call DDEST(cbuf,0) 0194 0195 C // Compute intermediate quantities 0196 0197 Coscl = cos(Colat0) 0198 Tancl = tan(Colat0) 0199 Tancl2 = tan(Colat0/2.) 0200 Mxtheta = Pi*Coscl 0201 0202 write(cbuf,'('' nvxini: Coscl, Tancl'', 2F7.4)') 0203 * Coscl, Tancl 0204 call DDEST(cbuf,0)
0205 write(cbuf,'('' nvxini: Tancl2, Mxtheta '', 2F7.4)') 0206 * tancl2, Mxtheta 0207 call DDEST(cbuf,0) 0208 0209 Latlon = .TRUE. 0210 0211 0212 else ! option=1 but type not 'TANC' 0213 0214 call DDEST('nvxini(tanc) parameter type bad',0) 0215 Init = .FALSE. 0216 NVXINI = -4 0217 return 0218 0219 end if 0220 0221 else if ( option .eq. 2) then 0222 0223 call movwc(param(1),outcoord) 0224 if( outcoord.eq.'LL' ) then 0225 Latlon = .TRUE. 0226 else if( outcoord.eq.'XYZ') then 0227 Latlon = .FALSE. 0228 else 0229 call DDEST('option=2 coord '//outcoord//' not supported',0) 0230 Init = .FALSE. 0231 NVXINI = -5 0232 end if 0233 0234 else ! option not 1 or 2 0235 0236 call DDEST('nvxini(tanc) unrecognized output option ',option) 0237 NVXINI = -4 0238 return 0239 0240 end if 0241 0242 NVXINI = 0 0243 Init = .TRUE. 0244 0245 return 0246 end 0247 0248 0249 0250 0251 *$ Name: 0252 *$ nvxsae - Compute earth coordinates from image coordinates 0253 *$ 0254 *$ Interface: 0255 *$ integer function 0256 *$ nvxsae( real lin, real ele, real dummy, 0257 *$ real e1, real e2, real e3 ) 0258 *$ 0259 *$ Input: 0260 *$ lin - image line 0261 *$ ele - image element 0262 *$ dummy - (unused) 0263 *$ 0264 *$ Input and Output: 0265 *$ none 0266 *$ 0267 *$ Output: 0268 *$ e1 - latitude or x 0269 *$ e2 - longitude or y 0270 *$ e3 - height or z 0271 *$ 0272 *$ Return values: 0273 *$ 0 - success 0274 *$ -1 - input data physically valid but not navigable 0275 *$ given the specified projection 0276 *$ -6 - module not initialized 0277 *$ 0278 *$ Remarks: 0279 *$ The navigation module must first be initialized with 0280 *$ a call to nvxini(). The output form (lat,lon) or (x,y,z) 0281 *$ depends on the last call to nvxini() with option 2. 0282 *$ 0283 *$ Categories: 0284 *$ navigation 0285 0286 INTEGER FUNCTION NVXSAE( lin, ele, dummy, e1, e2, e3 ) 0287 0288 implicit NONE 0289 0290 C // Interface variables (formal arguments) 0291 0292 real lin ! image line to navigate 0293 real ele ! image element to navigate 0294 real dummy ! (unused argument) 0295 real e1 ! Earth coordinate 1 0296 real e2 ! Earth coordinate 2 0297 real e3 ! Earth coordinate 3 0298 0299 C // Local variables 0300 0301 real lat ! latitude (McIDAS convention) 0302 real lon ! longitude (McIDAS convention) 0303 real hgt ! height 0304 real dx ! zonal displacement from pole 0305 C ! on projection surface 0306 real dy ! meridional displacement from pole 0307 real radius ! distance from pole on projection 0308 real theta ! angle from standard longitude on 0309 C ! projection surface 0310 real colat ! colatitude of navigated point 0311 0312 0313 C /////////////////////////////////////////////////////// 0314 C Common block variables and declaration. 0315 0316 C ALL CODE BETWEEN THE '//////' SEPARATORS MUST BE 0317 C DUPLICATED EXACTLY IN EACH NAVIGATION ROUTINE 0318 0319 C (A more maintenance-safe version would use ENTRY points 0320 C rather than separate functions for the navigation API 0321 C but entry points cannot be processed by 'convdlm.') 0322 0323 C // Common block contents: projection parameters 0324 0325 real Lin0 ! image line of pole 0326 real Ele0 ! image element of pole 0327 real Scale ! km per unit image coordinate 0328 C ! (pixel) 0329 real Lon0 ! standard longitude 0330 real Colat0 ! standard colatitude 0331 0332 C // Common block contents: pre-computed intermediate values 0333 0334 real Coscl ! cosine(Colat0) 0335 real Tancl ! tangent(Colat0) 0336 real Tancl2 ! tangent(Colat0/2) 0337 real Mxtheta ! limit of angle from std. lon 0338 C ! on projection surface 0339 0340 C // Common block contents: constants 0341 0342 real D2R ! degrees to radians factor 0343 real Pi 0344 real Badreal ! returned when navigation 0345 C ! cannot be done 0346 real Erad ! Earth radius 0347 logical Init ! initialized flag 0348 logical Latlon ! .TRUE. for lat/lon I/O 0349 0350 COMMON/TANC/ Lin0, Ele0, Scale, Lon0, Colat0, 0351 & Coscl, Tancl, Tancl2, Mxtheta, 0352 & D2R, Pi, Badreal, Erad, Init, Latlon 0353 0354 C End of common block variables and declaration. 0355 C /////////////////////////////////////////////////// 0356 0357 e1 = Badreal 0358 e2 = Badreal 0359 e3 = Badreal 0360 0361 0362 C // verify initialized module 0363 0364 if(.not.Init) then 0365 NVXSAE = -6 0366 return 0367 end if 0368 0369 0370 C // Compute radius and bearing from pole 0371 0372 dx = Scale*(lin-Lin0) 0373 dy = Scale*(ele-Ele0) 0374 0375 radius = sqrt(dx*dx+dy*dy) 0376 theta = atan2(dy,dx) 0377 0378 0379 C // Apply range checking on theta to determine if point is navigable 0380 0381 if ( theta.le.-Mxtheta .or. theta.gt.Mxtheta ) then 0382 NVXSAE = -1 0383 return 0384 end if 0385 0386 C // Forward navigation: compute longitude and colatitude 0387 C // from radius and theta 0388 0389 lon = Lon0 + theta/Coscl 0390 if(lon.le.-Pi) lon = lon + 2.d0*Pi 0391 if(lon.gt. Pi) lon = lon - 2.d0*Pi 0392 0393 colat = 2.*atan( Tancl2 * (radius/(Erad*Tancl))**(1./Coscl)) 0394 0395 C // Rescale to McIDAS convention (degrees, West positive). 0396 C // Apply conversion to Cartesian coordinates if 'XYZ' set 0397 C // as output form. Set return code for success. 0398 0399 lon = -lon/D2R 0400 lat = 90. - colat/D2R 0401 hgt = 0. 0402 0403 if(.not.Latlon) then 0404 call nllxyz(lat,lon,e1,e2,e3) 0405 else 0406 e1 = lat 0407 e2 = lon 0408 e3 = 0. 0409 end if 0410 0411 NVXSAE = 0 0412 0413 return 0414 end
0415 *$ Name: 0416 *$ nvxeas - Compute image coordinates from earth coordinates 0417 *$ 0418 *$ Interface: 0419 *$ integer function 0420 *$ nvxeas( real e1, real e2, real e3, 0421 *$ real lin, real ele, real dummy) 0422 *$ 0423 *$ Input: 0424 *$ e1 - latitude or x 0425 *$ e2 - longitude or y 0426 *$ e3 - height or z 0427 *$ 0428 *$ Input and Output: 0429 *$ none 0430 *$ 0431 *$ Output: 0432 *$ lin - image line 0433 *$ ele - image element 0434 *$ dummy - (unused) 0435 *$ 0436 *$ Return values: 0437 *$ 0 - success 0438 *$ -1 - input data physically valid but not navigable 0439 *$ given the specified projection 0440 *$ -2 - input data exceed physical limits 0441 *$ -6 - module not initialized 0442 *$ 0443 *$ Remarks: 0444 *$ The navigation module must first be initialized with 0445 *$ a call to nvxini(). The input form (lat,lon) or (x,y,z) 0446 *$ depends on the last call to nvxini() with option 2. 0447 *$ Input longitude may be in the range -360 to +360; 0448 *$ values outside this range will not be de-navigated. 0449 *$ Height (hgt) is ignored. 0450 *$ 0451 *$ Categories: 0452 *$ navigation 0453 0454 INTEGER FUNCTION NVXEAS( e1, e2, e3, lin, ele, dummy) 0455 0456 implicit NONE 0457 0458 C // Interface variables (formal arguments) 0459 0460 real e1 ! Earth coordinate 1 0461 real e2 ! Earth coordinate 2 0462 real e3 ! Earth coordinate 3 0463 real lin ! image line to navigate 0464 real ele ! image element to navigate 0465 real dummy ! (unused argument) 0466 0467 C // Local variables 0468 0469 real lat ! latitude (McIDAS convention) 0470 real lon ! longitude (McIDAS convention) 0471 real hgt ! height 0472 real in_lon ! input longitude (radians, 0473 C ! East positive) 0474 real colat ! colatitude 0475 real radius ! distance from pole on projection 0476 real theta ! angle from standard longitude on 0477 C ! projection surface 0478 0479 C //////////////////////////////////////////////////////// 0480 C Common block variables and declaration. 0481 0482 C ALL CODE BETWEEN THE '//////' SEPARATORS MUST BE 0483 C DUPLICATED EXACTLY IN EACH NAVIGATION ROUTINE 0484 0485 C (A more maintenance-safe version would use ENTRY points 0486 C rather than separate functions for the navigation API 0487 C but entry points cannot be processed by 'convdlm.') 0488 0489 C // Common block contents: projection parameters 0490 0491 real Lin0 ! image line of pole 0492 real Ele0 ! image element of pole 0493 real Scale ! km per unit image coordinate 0494 C ! (pixel) 0495 real Lon0 ! standard longitude 0496 real Colat0 ! standard colatitude 0497 0498 C // Common block contents: precomputed intermediate values 0499 0500 real Coscl ! cosine(Colat0) 0501 real Tancl ! tangent(Colat0) 0502 real Tancl2 ! tangent(Colat0/2) 0503 real Mxtheta ! limit of angle from std. lon 0504 C ! on projection surface 0505 0506 C // Common block contents: constants 0507 0508 real D2R ! degrees to radians factor 0509 real Pi 0510 real Badreal ! returned when navigation 0511 C ! cannot be done 0512 real Erad ! Earth radius 0513 logical Init ! initialized flag 0514 logical Latlon ! .TRUE. for lat/lon I/O 0515 0516 0517 COMMON/TANC/ Lin0, Ele0, Scale, Lon0, Colat0, 0518 & Coscl, Tancl, Tancl2, Mxtheta, 0519 & D2R, Pi, Badreal, Erad, Init, Latlon 0520 0521 C End of common block variables and declaration. 0522 C ///////////////////////////////////////////////////// 0523 0524 lin = Badreal 0525 ele = Badreal 0526 dummy = Badreal 0527 0528 C // verify that module is initialized 0529 0530 if(.not.init) then 0531 NVXEAS = -6 0532 return 0533 end if 0534 0535 C // Preprocess input values. If mode is 'XYZ' first convert 0536 C // from Cartesian to lat/lon. If mode is 'LL' just transcribe 0537 C // from arguments. 0538 0539 if(Latlon) then 0540 lat = e1 0541 lon = e2 0542 hgt = e3 0543 else 0544 call nxyzll( e1, e2, e3, lat, lon) 0545 hgt = 0. 0546 end if 0547 0548 C // check that input values are physically possible and 0549 C // then convert to radians and East positive 0550 0551 if ( lat.lt.-90. .or. lat.gt.90. ) then 0552 NVXEAS = -2 0553 return 0554 end if
0555 if( lon.le.-360..or.lon.gt.360.) then 0556 NVXEAS = -2 0557 return 0558 end if 0559 0560 if( lat.eq.-90. .or. lat.eq.90. ) then 0561 NVXEAS = -1 0562 return 0563 end if 0564 0565 colat = Pi/2. - D2R*lat 0566 in_lon = -D2R*lon 0567 0568 C // map longitude into range -Pi to Pi 0569 0570 if(in_lon.le.-Pi) in_lon = in_lon + 2.*Pi 0571 if(in_lon.gt. Pi) in_lon = in_lon - 2.*Pi 0572 0573 0574 C // Now trap South Pole. Though a physically possible latitude, 0575 C // tan(colat/2) -> infinity there so it is not navigable 0576 0577 if ( colat.eq.Pi ) then 0578 NVXEAS = -1 0579 return 0580 end if 0581 0582 0583 C // Compute radius and theta of point on projection surface. 0584 C // Theta is tricky; you have to compute offset relative 0585 C // to standard longitude, force that into -pi to +pi range, 0586 C // and THEN scale by cos(Colat0) 0587 0588 radius = Erad * Tancl *( tan(colat/2.)/Tancl2 ) ** Coscl 0589 theta = in_lon-Lon0 0590 0591 if(theta.le.-Pi) theta = theta + 2.*Pi 0592 if(theta.gt. Pi) theta = theta - 2.*Pi 0593 0594 theta = Coscl * theta 0595 0596 0597 C // Compute line and element 0598 0599 lin = Lin0 + radius*cos(theta)/Scale 0600 ele = Ele0 + radius*sin(theta)/Scale 0601 dummy = 0. 0602 0603 NVXEAS = 0 0604 0605 return 0606 end 0607 0608 0609 *$ Name: 0610 *$ nvxopt - Perform supplemental navigation operations 0611 *$ 0612 *$ Interface: 0613 *$ integer function 0614 *$ nvxopt(integer option, real xin(*), 0615 *$ real xout(*) ) 0616 *$ Input: 0617 *$ option - 'SCAL' compute projection scale 0618 *$ xin(1) - latitude 0619 *$ 0620 *$ Input and Output: 0621 *$ none 0622 *$ 0623 *$ Output: 0624 *$ xout(1) - km per pixel at given latitude
0625 *$ Return values: 0626 *$ 0 - success 0627 *$ -1 - input latitude physically valid, but projection 0628 *$ undefined or scale infinite there 0629 *$ -2 - input latitude exceeds physical limits 0630 *$ -5 - unrecognized option 0631 *$ -6 - module not initialized 0632 *$ 0633 *$ Remarks: 0634 *$ The navigation module must first be initialized by 0635 *$ a call to nvxini(). Latitude is in degrees, north positive, 0636 *$ and must lie between -90. and +90. 0637 *$ 0638 *$ Categories: 0639 *$ navigation 0640 0641 INTEGER FUNCTION NVXOPT( option, xin, xout) 0642 0643 implicit NONE 0644 0645 C // Interface variables (formal arguments) 0646 0647 integer option ! special service name (character 0648 C ! stored as integer) 0649 real xin(*) ! input vector 0650 real xout(*) ! output vector 0651 0652 C // Local variables 0653 0654 character*4 copt ! special service (character form) 0655 real colat ! input colatitude 0656 0657 C //////////////////////////////////////////////////// 0658 C Common block variables and declaration. 0659 0660 C ALL CODE BETWEEN THE '//////' SEPARATORS MUST BE 0661 C DUPLICATED EXACTLY IN EACH NAVIGATION ROUTINE 0662 0663 C (A more maintenance-safe version would use ENTRY points 0664 C rather than separate functions for the navigation API 0665 C but entry points cannot be processed by 'convdlm.') 0666 0667 C // Common block contents: projection parameters 0668 0669 real Lin0 ! image line of pole 0670 real Ele0 ! image element of pole 0671 real Scale ! km per unit image coordinate 0672 C ! (pixel) 0673 real Lon0 ! standard longitude 0674 real Colat0 ! standard colatitude 0675 0676 C // Common block contents: precomputed intermediate values 0677 0678 real Coscl ! cosine(Colat0) 0679 real Tancl ! tangent(Colat0) 0680 real Tancl2 ! tangent(Colat0/2) 0681 real Mxtheta ! limit of angle from std. lon 0682 C ! on projection surface 0683 0684 C // Common block contents: constants 0685 0686 real D2R ! degrees to radians factor 0687 real Pi 0688 real Badreal ! returned when navigation 0689 C ! cannot be done 0690 real Erad ! Earth radius 0691 logical Init ! initialized flag 0692 logical Latlon ! .TRUE. for lat/lon I/O 0693 0694 0695 COMMON/TANC/ Lin0, Ele0, Scale, Lon0, Colat0, 0696 & Coscl, Tancl, Tancl2, Mxtheta, 0697 & D2R, Pi, Badreal, Erad, Init, Latlon 0698 0699 C End of common block variables and declaration. 0700 C //////////////////////////////////////////////////// 0701 0702 0703 xout(1) = Badreal 0704 0705 C // verify initialized module 0706 0707 if(.not.init) then 0708 NVXOPT = -6 0709 return 0710 end if 0711 0712 C // Extract and interpret the option 0713 0714 call movwc(option,copt) 0715 0716 if(copt.eq.'SCAL') then 0717 0718 C // Compute colatitude and make sure it is 0719 C // physically possible and navigable 0720 0721 if ( xin(1).gt.90. .or. xin(1).lt.-90. ) then 0722 NVXOPT = -2 0723 return 0724 else if ( xin(1).eq.90. .or. xin(1).eq.-90. ) then 0725 NVXOPT = -1 0726 return 0727 end if 0728 0729 colat = Pi/2. - D2R*xin(1) 0730 0731 C // Now compute actual scale for this colatitude 0732 0733 xout(1) = scale 0734 * *(sin(Colat0)*(tan(colat/2.)/Tancl2)**Coscl)/sin(colat) 0735 0736 C else if(copt.eq.'????') 0737 C // Add code for additional options here 0738 0739 else 0740 NVXOPT = -5 0741 return 0742 end if 0743 0744 NVXOPT = 0 0745 0746 return 0747 end
Calibration is the process of converting data values sensed by an instrument to useful, physical quantities such as temperature, radiance and albedo. The McIDAS-X calibration subsystem is designed to:
Thus, you can define your own calibration modules that will allow McIDAS-X applications to view the data that you prescribe.
McIDAS-X has a defined Application Program Interface, including subroutine names, calling sequence and functionality, that all calibration modules must adhere to. This section provides a brief history of how the calibration API was developed and an overview of the applications interface.
In the mid-1980s, the McIDAS-X area-format image file was redesigned to accommodate the ever increasing number of remotely sensed data types. The redesign allowed programmers to add new datasets from a variety of platforms to McIDAS-X without changing the existing software or file format. The new area file format accommodated multibanded, multibyte data along with a variety of ancillary data. The redesign process also provided a general method of storing, navigating and calibrating these data. It defined an API that all navigation and calibration modules would adhere to, and a mechanism for accessing the appropriate module at application runtime. The previous section in this chapter described the navigation system; this section describes the calibration system.
|
For information about the McIDAS-X area format image file, see Chapter 6, Format of the Data Files. |
The calibration API is defined at several levels. For ADDE applications, the call to mcaget has a parameter to request data to be returned in a specific physical quantity. The calibration process occurs on the data server, and although calibration information is returned, it is usually not needed by the ADDE client application.
The next level down in the API is used by non-ADDE applications and by the ADDE data servers. A call to the araopt function, which sets the area options, will use a calibration module if the physical quantity specified in its UNIT option is different than that stored in the image file (word 53 of the area directory). Subsequent calls to redara , which reads the data, will use the calibration module if needed. This provides the application with calibrated data in a rather transparent and data-independent way.
Below araopt in the API is the subroutine kbprep , which provides the interface to the dynamic calibration modules. Most applications don't call the lower level API functions directly. However, you must understand this lowest API definition to incorporate new calibration modules.
The function names used to access the calibration are listed below:
These are not the names of the functions within the calibration module itself. A mapping is done in kbprep that allows applications to use multiple instances of different calibrations. kbprep builds the name of the module to load, based on the slot number and data type.
The numeric (1, 2 or 3) is the slot number, which allows simultaneous use of up to three different calibration modules. The data type, which is stored in word 52 of the area directory, is needed to construct the name of the module to call. For example, the data type for GOES-8 is GVAR. Thus, the name constructed, using slot number 1, is KB1GVAR, which uses the functions KB1INIGVAR, KB1CALGVAR, and KB1OPTGVAR.
All calibration modules have the same framework. They must conform to the McIDAS-X convention for functionality, names of functions, and types of arguments. This standardization allows applications to make use of these modules in a generic, yet powerful way. It is usually not necessary for applications to have any private knowledge of the data it's working with; the calibration interface provides a means to acquire certain aspects that are common to all.
When naming your calibration module, use this form: kbxtype.dlm, where type is the same four characters in word 52 of the area directory. For example, type could be GVAR, TIRO or AAA.
(1)
(2)
(3)
(4)
(5)
(6)
(7)
(8)
(9)
