[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: 20050926: perl/netCDF record mismatch
- Subject: Re: 20050926: perl/netCDF record mismatch
- Date: Mon, 26 Sep 2005 13:55:57 -0600 (MDT)
andres,
first the error you got is the result of writing the report to the netCDF
file. i tried to reproduce the error using the latest ua2nc and raob.cdl file
on data set 2005092612_upa.wmo, but it didn't produce any errors. i did
fix a similar error some time ago, but i don't remember the time period.
i'm going to attach a raob.cdl and ua2nc files for you to try and see if
you still get the error. after you install raob.cdl and ua2nc files, stop
the ldm and remove the files uaLog.log and ua.lst in the data directory or
you will get no data for one day. i'm getting ready to make a new release
so if you still have problems, i would like to know. also, there is a new
file naming convention, ie Upperair_20050926_0000.nc on a daily basis. if
you want the old naming convention, then add flag!! -n old. if you want
hourly files, then add -h. for old files and hourly, then:
% ua2nc -h -n old raob.cdl ....
robb...
Opening data/Upperair_20050926_0000.nc with ncid 3
Station: 91592, report type: TTAA, j = 12
Time: 2005092612
---------------------------------------------------------
Pres Ht Temp DewPt Dir Spd
mb m C C deg m/s
---------------------------------------------------------
1014 69 293.95 3 80 7.2016
1000 172 293.15 2.8 80 7.716
925 840 288.75 2.2 70 11.8312
850 1553 283.75 1.2 55 13.3744
700 3159 281.55 47 200 5.144
500 5850 265.45 47 285 10.288
400 7540 251.05 42 250 6.6872
300 9580 238.05 37 265 33.9504
250 10850 233.85 36 285 48.868
200 12330 222.05 99999 285 49.3824
150 14140 207.45 99999 290 43.724
100 16530 197.45 99999 290 21.6048
Maximum values
---------------------------------------------------------
226 285 51.44
Station: 91592, report type: TTBB, j = 9
Time: 2005092612
---------------------------------------------------------
Pres Ht Temp DewPt Dir Spd
mb m C C deg m/s
---------------------------------------------------------
1014.0 99999 293.95 3 99999 99999
794.0 99999 279.95 0.7 99999 99999
781.0 99999 278.95 5 99999 99999
758.0 99999 283.95 19 99999 99999
742.0 99999 282.55 24 99999 99999
713.0 99999 282.55 47 99999 99999
501.0 99999 265.65 47 99999 99999
384.0 99999 248.45 41 99999 99999
336.0 99999 239.85 14 99999 99999
On Mon, 26 Sep 2005, Unidata Support wrote:
>
> ------- Forwarded Message
>
> >To: address@hidden
> >From: =?ISO-8859-1?Q?Andr=E9s_Calder=F3n?= <address@hidden>
> >Subject: perl/netCDF record mismatch
> >Organization: ?
> >Keywords: 200509261616.j8QGGeYJ010181 netCDF decoders
>
> Hi,
>
> I am trying convert upperair data to netcdf, the command tested:
>
> $ ua2nc -v /usr/local/etc/raob.cdl <05092612_upa.wmo
>
> and the log is (an extract):
> ....
> Station: 91592, report type: TTAA, j =3D 12
> Time: 2005092612
> ---------------------------------------------------------
> Pres Ht Temp DewPt Dir Spd
> mb m C C deg m/s
> ---------------------------------------------------------
> 1014 69 293.95 3 80 7.2016
> 1000 172 293.15 2.8 80 7.716
> 925 840 288.75 2.2 70 11.8312
> 850 1553 283.75 1.2 55 13.3744
> 700 3159 281.55 47 200 5.144
> 500 5850 265.45 47 285 10.288
> 400 7540 251.05 42 250 6.6872
> 300 9580 238.05 37 265 33.9504
> 250 10850 233.85 36 285 48.868
> 200 12330 222.05 99999 285 49.3824
> 150 14140 207.45 99999 290 43.724
> 100 16530 197.45 99999 290 21.6048
>
> Maximum values
> ---------------------------------------------------------
> 226 285 51.44
>
> ncvarget: ncid 3: Edge+start exceeds dimension bound
> NetCDF::recget status =3D -1
> perl/netCDF record mismatch at /usr/local/bin/ua2nc line 355, <STDIN> chunk
> 7.
> Wmo =3D91592, Type =3DTTAA, Rpt: 2611, Trans: 2612, recput =3D-1
> .....
>
> What happens?
>
> What means? :
> ncvarget: ncid 3: Edge+start exceeds dimension bound
> NetCDF::recget status =3D -1
> perl/netCDF record mismatch at /usr/local/bin/ua2nc line 355, <STDIN> chunk
> 7.
>
>
>
> decoder 3.1.2
> debian sarge
> netcdf 3.5.0
> perl v5.8.7
>
> thanks,
>
> -- Andres Calderon
>
> --
> NOTE: All email exchanges with Unidata User Support are recorded in the
> Unidata inquiry tracking system and then made publicly available
> through the web. If you do not want to have your interactions made
> available in this way, you must let us know in each email you send to us.
>
> ------- End of Forwarded Message
>
===============================================================================
Robb Kambic Unidata Program Center
Software Engineer III Univ. Corp for Atmospheric Research
address@hidden WWW: http://www.unidata.ucar.edu/
===============================================================================
netcdf raob {
dimensions:
recNum = UNLIMITED ;
manLevel = 22 ;
sigTLevel = 150 ;
sigWLevel = 76 ;
mWndNum = 4 ;
mTropNum = 4 ;
staNameLen = 6 ;
variables:
long wmoStaNum (recNum) ;
wmoStaNum:long_name = "WMO Station Number" ;
wmoStaNum:reference = "Volume A of WMO publication 9" ;
wmoStaNum:_FillValue = 99999 ;
char staName (recNum, staNameLen) ;
staName:long_name = "Station Identifier" ;
float staLat (recNum) ;
staLat:long_name = "Station Latitude" ;
staLat:valid_range = 0.f, 90.f ;
staLat:_FillValue = 99999.f ;
staLat:units = "degrees_north" ;
float staLon (recNum) ;
staLon:long_name = "Station Longitude" ;
staLon:valid_range = -180.f, -50.f ;
staLon:_FillValue = 99999.f ;
staLon:units = "degrees_east" ;
float staElev (recNum) ;
staElev:long_name = "Station Elevation" ;
staElev:valid_range = -100.f, 3500.f ;
staElev:_FillValue = 99999.f ;
staElev:units = "meter" ;
double synTime (recNum) ;
synTime:long_name = "Synoptic Time" ;
synTime:units = "seconds since (1970-1-1 00:00:0.0)" ;
long numMand (recNum) ;
numMand:long_name = "Number of Mandatory Levels" ;
numMand:valid_range = 0, 22 ;
long numSigT(recNum) ;
numSigT:long_name = "Number of Significant Levels wrt T" ;
numSigT:valid_range = 0, 150 ;
long numSigW(recNum) ;
numSigW:long_name = "Number of Significant Levels wrt W" ;
numSigW:valid_range = 0, 76 ;
long numMwnd(recNum) ;
numMwnd:long_name = "Number of Maximum Wind Levels" ;
numMwnd:valid_range = 0, 4 ;
long numTrop(recNum) ;
numTrop:long_name = "Number of Tropopause Levels" ;
numTrop:valid_range = 0, 4 ;
double relTime(recNum) ;
relTime:long_name = "Sounding Release Time" ;
relTime:_FillValue = -1.0e+38d;
relTime:units = "seconds since (1970-1-1 00:00:0.0)" ;
long sondTyp(recNum) ;
sondTyp:long_name = "Instrument Type" ;
sondTyp:reference = "Federal Meteorological Handbook No. 4" ;
sondTyp:_FillValue = 99999 ;
float prMan(recNum, manLevel) ;
prMan:long_name = "Pressure - Mandatory level" ;
prMan:valid_range = 1.f, 1500.f ;
prMan:_FillValue = 99999.f ;
prMan:units = "hectopascal" ;
float htMan(recNum, manLevel) ;
htMan:long_name = "Geopotential - Mandatory level" ;
htMan:valid_range = -250.f, 60000.f ;
htMan:_FillValue = 99999.f ;
htMan:units = "meter" ;
float tpMan(recNum, manLevel) ;
tpMan:long_name = "Temperature - Mandatory level" ;
tpMan:valid_range = 173.f, 373.f ;
tpMan:_FillValue = 99999.f ;
tpMan:units = "kelvin" ;
float tdMan(recNum, manLevel) ;
tdMan:long_name = "Dew Point Depression - Mandatory level" ;
tdMan:valid_range = 0.f, 60.f ;
tdMan:_FillValue = 99999.f ;
tdMan:units = "kelvin" ;
float wdMan(recNum, manLevel) ;
wdMan:long_name = "Wind Direction - Mandatory level" ;
wdMan:valid_range = 0.f, 360.f ;
wdMan:_FillValue = 99999.f ;
wdMan:units = "degree_true" ;
float wsMan(recNum, manLevel) ;
wsMan:long_name = "Wind Speed - Mandatory level" ;
wsMan:valid_range = 0.f, 300.f ;
wsMan:_FillValue = 99999.f ;
wsMan:units = "meter/sec" ;
float prSigT(recNum, sigTLevel) ;
prSigT:long_name = "Pressure - Significant level wrt T" ;
prSigT:valid_range = 1.f, 1500.f ;
prSigT:_FillValue = 99999.f ;
prSigT:units = "hectopascal" ;
float tpSigT(recNum, sigTLevel) ;
tpSigT:long_name = "Temperature - Significant level wrt T" ;
tpSigT:valid_range = 173.f, 373.f ;
tpSigT:_FillValue = 99999.f ;
tpSigT:units = "kelvin" ;
float tdSigT(recNum, sigTLevel) ;
tdSigT:long_name = "Dew Point Depression - Significant level
wrt T" ;
tdSigT:valid_range = 0.f, 60.f ;
tdSigT:_FillValue = 99999.f ;
tdSigT:units = "kelvin" ;
float htSigW(recNum, sigWLevel) ;
htSigW:long_name = "Geopotential - Significant level wrt W" ;
htSigW:valid_range = -250.f, 60000.f ;
htSigW:_FillValue = 99999.f ;
htSigW:units = "meter" ;
float wdSigW(recNum, sigWLevel) ;
wdSigW:long_name = "Wind Direction - Significant level wrt W" ;
wdSigW:valid_range = 0.f, 360.f ;
wdSigW:_FillValue = 99999.f ;
wdSigW:units = "degree_true" ;
float wsSigW(recNum, sigWLevel) ;
wsSigW:long_name = "Wind Speed - Significant level wrt W" ;
wsSigW:valid_range = 0.f, 300.f ;
wsSigW:_FillValue = 99999.f ;
wsSigW:units = "meter/sec" ;
float prTrop(recNum, mTropNum) ;
prTrop:long_name = "Pressure - Tropopause level" ;
prTrop:valid_range = 1.f, 1500.f ;
prTrop:_FillValue = 99999.f ;
prTrop:units = "hectopascal" ;
float tpTrop(recNum, mTropNum) ;
tpTrop:long_name = "Temperature - Tropopause level" ;
tpTrop:valid_range = 173.f, 373.f ;
tpTrop:_FillValue = 99999.f ;
tpTrop:units = "kelvin" ;
float tdTrop(recNum, mTropNum) ;
tdTrop:long_name = "Dew Point Depression - Tropopause level" ;
tdTrop:valid_range = 0.f, 60.f ;
tdTrop:_FillValue = 99999.f ;
tdTrop:units = "kelvin" ;
float wdTrop(recNum, mTropNum) ;
wdTrop:long_name = "Wind Direction - Tropopause level" ;
wdTrop:valid_range = 0.f, 360.f ;
wdTrop:_FillValue = 99999.f ;
wdTrop:units = "degree_true" ;
float wsTrop(recNum, mTropNum) ;
wsTrop:long_name = "Wind Speed - Tropopause level" ;
wsTrop:valid_range = 0.f, 300.f ;
wsTrop:_FillValue = 99999.f ;
wsTrop:units = "meter/sec" ;
float prMaxW(recNum, mWndNum) ;
prMaxW:long_name = "Pressure - Maximum wind level" ;
prMaxW:valid_range = 1.f, 1500.f ;
prMaxW:_FillValue = 99999.f ;
prMaxW:units = "hectopascal" ;
float wdMaxW(recNum, mWndNum) ;
wdMaxW:long_name = "Wind Direction - Maximum wind level" ;
wdMaxW:valid_range = 0.f, 360.f ;
wdMaxW:_FillValue = 99999.f ;
wdMaxW:units = "degree_true" ;
float wsMaxW(recNum, mWndNum) ;
wsMaxW:long_name = "Wind Speed - Maximum wind level" ;
wsMaxW:valid_range = 0.f, 300.f ;
wsMaxW:_FillValue = 99999.f ;
wsMaxW:units = "meter/sec" ;
// global attributes:
:comment0 = "First mandatory level is surface level" ;
:version = "Forecast Systems Lab 1.4" ;
}
#!/usr/bin/perl
#
# Decodes Upperair (WMO TEMP) reports according to:
# WMO Manual on codes volume 1 No. 306
# FM 35-IX Ext. TEMP - Upper-level pressure, temperture, humidity, and
# wind reports from a fixed land station.
#
# usage with input from STDIN: ua2nc [-v] [-h] [-n (old|new)] cdlfile
[datadir] [yyyymm]
#
# -v flag prints out the reports in the uaLog.log file, very verbose
#
use NetCDF ;
use Time::Local ;
#no encoding;
# process command line switches
$upperair = "new";
while( $_ = $ARGV[0], /^-/ ) {
shift ;
last if /^--$/ ;
/^(-v)/ && $verbose++ ;
/^(-h)/ && $hourly++; #create hourly files instead of daily
/^(-n)/ && ( $upperair = shift );
}
# new netCDF file naming covention is default
if( $upperair eq "new" ) {
if( $hourly ) {
$upperair = "00";
} else {
$upperair = "_0000";
}
$upperairPrefix = "Upperair_";
} else {
$upperair = "_upperair";
$upperairPrefix = "";
}
# process input parameters
if( $#ARGV == 0 ) {
$cdlfile = $ARGV[ 0 ] ;
} elsif( $#ARGV == 1 ) {
$cdlfile = $ARGV[ 0 ] ;
if( $ARGV[ 1 ] =~ /^\d/ ) {
$yyyymm = $ARGV[ 1 ] ;
} else {
$datadir = $ARGV[ 1 ] ;
}
} elsif( $#ARGV == 2 ) {
$cdlfile = $ARGV[ 0 ] ;
$datadir = $ARGV[ 1 ] ;
$yyyymm = $ARGV[ 2 ] ;
} else {
die "Usage with input on STDIN: ua2nc [-v] [-h] [-n (old|new)] cdlfile
[datadir] [yyyymm]\n" ;
}
# check for cdl and netCDF ncgen
die "Missing cdlfile parameter: $!\n" unless -e $cdlfile ;
if( -e "util/ncgen" ) {
$ncgen = "util/ncgen" ;
} elsif( -e "/usr/local/ldm/util/ncgen" ) {
$ncgen = "/usr/local/ldm/util/ncgen" ;
} elsif( -e "/upc/netcdf/bin/ncgen" ) {
$ncgen = "/upc/netcdf/bin/ncgen" ;
} elsif( -e "./ncgen" ) {
$ncgen = "./ncgen" ;
} else {
open( NCGEN, "which ncgen |" ) ;
$ncgen = <NCGEN> ;
close( NCGEN ) ;
if( $ncgen =~ /no ncgen/ ) {
die "Can't find NetCDF utility 'ncgen' in PATH, util/ncgen
/usr/local/ldm/util/ncgen, /upc/netcdf/bin/ncgen, or ./ncgen : $!\n" ;
} else {
$ncgen = "ncgen" ;
}
}
# the data directory
$datadir = "." if( ! $datadir ) ;
system( "mkdir -p $datadir" ) if( ! -e $datadir ) ;
if( -e "$datadir/uaLog.log.2" ) {
`rm -f $datadir/uaLog.log.3`;
`mv $datadir/uaLog.log.2 $datadir/uaLog.log.3`;
}
if( -e "$datadir/uaLog.log.1" ) {
`mv $datadir/uaLog.log.1 $datadir/uaLog.log.2`;
}
if( -e "$datadir/uaLog.log" ) {
`mv $datadir/uaLog.log $datadir/uaLog.log.1`;
}
# redirect STDOUT and STDERR
open( STDOUT, ">$datadir/uaLog.log" ) ||
die "could not open $datadir/uaLog.log: $!\n" ;
open( STDERR, ">&STDOUT" ) ||
die "could not dup stdout: $!\n" ;
select( STDERR ) ; $| = 1 ;
select( STDOUT ) ; $| = 1 ;
# year and month
if( ! $yyyymm ) {
$theyear = (gmtime())[ 5 ] ;
$theyear = ( $theyear < 100 ? $theyear : $theyear - 100 ) ;
$thedecade = sprintf( "%02d", $theyear ) ;
$theyear = "20" . sprintf( "%02d", $theyear ) ;
$themonth = (gmtime())[ 4 ] ;
$themonth++ ;
$themonth = sprintf( "%02d", $themonth ) ;
$yyyymm = "$theyear$themonth" ;
} else {
die "yyyymm must be 6 in length: $!\n" if( length( $yyyymm ) != 6 ) ;
$theyear = substr( $yyyymm, 0, 4 ) ;
$themonth = substr( $yyyymm, 4 ) ;
$thedecade = substr( $yyyymm, 2, 2 ) ;
}
%bin = ( "23", "00", "00", "00", "01", "00", "02", "03", "03", "03",
"04", "03", "05", "06", "06", "06", "07", "06", "08", "09",
"09", "09", "10", "09", "11", "12", "12", "12", "13", "12",
"14", "15", "15", "15", "16", "15", "17", "18", "18", "18",
"19", "18", "20", "21", "21", "21", "22", "21" ) ;
# set error handling to verbose only
$status = NetCDF::opts( VERBOSE ) ;
# set interrupt handler
$SIG{ 'INT' } = 'atexit' ;
$SIG{ 'KILL' } = 'atexit' ;
$SIG{ 'TERM' } = 'atexit' ;
$SIG{ 'QUIT' } = 'atexit' ;
# set defaults
%mand_ttaa = ("99", 99999, "00", 1000, "92", 925, "85", 850, "70", 700,
"50", 500, "40", 400, "30", 300, "25", 250, "20", 200,
"15", 150, "10", 100 ) ;
%mand_ttcc = ("70", 1, "50", 1, "30", 1, "20", 1, "10", 1, "07", 1, "05", 1,
"03", 1, "02", 1, "01", 1 ) ;
%mand_pres = ( "1000", 1, "925", 2, "850", 3, "700", 4, "500", 5, "400", 6,
"300", 7, "250", 8, "200", 9, "150", 10, "100", 11 , "70", 12,
"50", 13, "30", 14, "20", 15, "10", 16, "07", 17, "05", 18, "03", 19,
"02", 20, "01", 21 ) ;
%TTAA_top_wind_level = ( "/", 99999, "0", 1000, "9", 925, "8", 850, "7", 700,
"5", 500, "4", 400, "3", 300, "2", 200, "1", 100 ) ;
%TTCC_top_wind_level = ( "/", 100, "0", 0, "7", 70, "5", 50, "3", 30,
"2", 20, "1", 0) ;
# defaults
$manLevel = 22 ;
$sigTLevel = 65 ;
$sigWLevel = 50 ;
$mWndNum = 4 ;
$mTropNum = 4 ;
$staNameLen = 6 ;
# open cdl and get dimensions for variables
open( CDL, "$cdlfile" ) || die "could not open $cdlfile: $!\n" ;
$i = 0 ;
while( <CDL> ) {
if( s#^\s*manLevel\s*=\s*(\d{1,5})## ) {
$manLevel = $1 ;
} elsif( s#^\s*sigTLevel\s*=\s*(\d{1,5})## ) {
$sigTLevel = $1 ;
} elsif( s#^\s*sigWLevel\s*=\s*(\d{1,5})## ) {
$sigWLevel = $1 ;
} elsif( s#^\s*mWndNum\s*=\s*(\d{1,5})## ) {
$mWndNum = $1 ;
} elsif( s#^\s*mTropNum\s*=\s*(\d{1,5})## ) {
$mTropNum = $1 ;
} elsif( s#^\s*staNameLen\s*=\s*(\d{1,5})## ) {
$staNameLen = $1 ;
} elsif( s#^\s*variables## ) {
last ;
}
}
close CDL ;
# initialize record buffer and variables
makedataref();
# read in station data
if( -e "etc/snworld.tbl" ) {
$sfile = "etc/snworld.tbl" ;
} elsif( -e "./snworld.tbl" ) {
$sfile = "./snworld.tbl" ;
} else {
die "Can't find snworld.tbl station file.: $!\n" ;
}
open( STATION, "$sfile" ) || die "could not open $sfile: $!\n" ;
while( <STATION> ) {
s#^(\w{3,6})?\s+(\d{4,5}).{40}## ;
$id = $1 ;
$wmo_id = $2 ;
$wmo_id = "0" . $wmo_id if( length( $wmo_id ) == 4 ) ;
( $lat, $lon, $elev ) = split ;
$lat = sprintf( "%7.2f", $lat / 100 ) ;
$lon = sprintf( "%7.2f", $lon / 100) ;
# set these vars ( $lat, $lon, $elev, $id )
$STATIONS{ "$wmo_id" } = "$lat $lon $elev $id" ;
}
close STATION ;
# read in list of already processed reports if it exists
# open ua.lst, list of reports processed in the last 24 hours.
if( -e "$datadir/ua.lst" ) {
open( LST, "$datadir/ua.lst" ) ||
die "could not open $datadir/ua.lst: $!\n" ;
while( <LST> ) {
( $stn, $yyyymmddhh, $record, @R ) = split ;
$rpt_hash{ "$stn $yyyymmddhh" } = "$record @R" ;
}
close LST ;
}
# Now begin parsing file and decoding observations breaking on cntrl C
$/ = "\cC" ;
# set select processing here from STDIN
START:
while( 1 ) {
open( STDIN, '-' ) ;
vec( $rin,fileno( STDIN ),1 ) = 1 ;
$timeout = 1200 ; # 20 minutes
$nfound = select( $rout = $rin, undef, undef, $timeout ) ;
# timed out
if( ! $nfound ) {
print "Shut down, time out $timeout seconds on select\n" ;
atexit() ;
}
atexit( "eof" ) if( eof( STDIN )) ;
# Process each line of ua bulletins, header first
$_ = <STDIN> ;
next unless( /TT|PP/ ) ;
s#\cC## ;
s#\cM##g ;
s#\cA\n## ;
s#\c^##g ;
s#^\n##g ;
s#\d\d\d \n## ;
s#\w{4}\d{1,2} \w{4} (\d{2})(\d{2})(\d{2}).*\n## ;
$bday = $1 ;
$bhour = $2 ;
$bhour = "23" if( $bhour eq "24" ) ;
$bmin = $3 ;
next unless ( $bday && defined( $bhour ) && defined( $bmin ) ) ;
next if( $bmin > 59 || $bhour > 23 || $bday > 31 ) ;
# check for valid transmission times against current time
$cday = (gmtime())[ 3 ] ;
$chour = (gmtime())[ 2 ] ;
# skip bulletins over 24 hours old or in the future
if( $bday == $cday ) {
next if( $bhour > $chour ) ;
} else { # $bday != $cday, skip over day old reports
next if( $bday < $cday -1 ) ;
if( $bday > $cday ) {
next if( $cday != 1 || $bday < 28) ;
}
next if( $bhour < $chour ) ;
}
# Separate bulletins into reports
if( /=\n/ ) {
s#=\s+\n#=\n#g ;
} else {
s#\n# #g ;
}
@reports = split( /=\n/ ) ;
for( @reports ) { # Process each report in the bulletin
undef( $stn_name ) ;
next if( /NIL/ ) ;
s#\n# #g ;
s#^\d{5}.*((TT)(AA|BB|CC|DD)|(PP)(BB|DD))#$1# ;
next unless (
s# ?((TT)(AA|BB|CC|DD)|(PP)(BB|DD)) {1,2}(\d\d)(\d\d)([\d\/]) (\d{5}) ## ) ;
#print "$&\n" ;
$type = $1 ;
$rday = $6 ;
$rhour = $7 ;
$indicator = $8 ;
$stn = $9 ;
if( $rday =~ /^[5678]/ ) {
$rday -= 50 ;
$knots = 1 ;
} else {
$knots = 0 ;
}
# check for valid times
next unless ($rday && defined( $rhour )) ;
next if( $rhour > 23 || $rday > 31 ) ;
# skip reports over 24 hours old
$tmpyyyymm = $yyyymm ;
if( $rday == $cday ) {
next if( $rhour > $chour ) ;
} else { # $rday != $cday, skip over day old reports
next if( $rday < ( $cday -1 ) ) ;
if( $rday > $cday ) {
next if( $cday != 1 || $rday < 28) ;
# cday = 1, reset month and year
$tmpmonth = sprintf( "%02d", $themonth -1 ) ;
if( $tmpmonth == 0 ) {
$tmpmonth = "12" ;
$tmpyear =
sprintf( "%04d", $theyear -1 );
} else {
$tmpyear = $theyear ;
}
$tmpyyyymm = $tmpyear . $tmpmonth ;
}
next if( $rhour < $chour ) ;
}
$yyyymmddhh = $tmpyyyymm . sprintf( "%02d", $rday ) .
$bin{ "$rhour" } ;
if( defined( $STATIONS{ "$stn" } )) {
# station's lat, lon, elev, id.
( $lat, $lon, $elev, $id ) =
split( ' ', $STATIONS{ "$stn" } ) ;
$id = $AS6 unless( defined( $id )) ;
} else {
print "WMO Id $stn not in station table\n" ;
next ;
}
if( defined( $rpt_hash{ "$stn $yyyymmddhh" } )) {
( $record, @R ) = split( ' ',
$rpt_hash{ "$stn $yyyymmddhh" } ) ;
$rpt_hash{ "$stn $yyyymmddhh" } .= " $type"
if( index( "@R", $type ) == -1 ) ;
} else {
$record = -1 ;
}
# clean data
s#\s+# #g ;
$_ .= ' ' ;
$j = 0 ;
# process reports
SWITCH: {
ttaa( $_ ), last SWITCH if( $type eq "TTAA" ) ;
ttbb( $_ ), last SWITCH if( $type eq "TTBB" ) ;
ttcc( $_ ), last SWITCH if( $type eq "TTCC" ) ;
ttdd( $_ ), last SWITCH if( $type eq "TTDD" ) ;
ppbbdd( $_ ), last SWITCH
if( $type eq "PPBB" || $type eq "PPDD" ) ;
}
# no data extracted
next unless $j ;
# covert knots to meters/second if necessary
if( $knots && $#wspd != -1 ) {
for( $i = 0; $i < $j; $i++ ) {
$wspd[ $i ] *= 0.5144 if( $wspd[ $i ] != $F ) ;
}
}
# we have a legal report, open a Netcdf file
$status = doNet() ;
# print output
printvars() if( $verbose ) ;
# output the NetCDF data here
makedataref();
# record already exists for $stn
if( $record >= 0 ) {
$status = NetCDF::recget( $ncid, $record, \@dataref ) ;
print "NetCDF::recget status = $status\n"
if( $status ) ;
} else { # initialize new record
$record = $recnum ;
$rpt_hash{ "$stn $yyyymmddhh" } = "$record $type" ;
$recnum++ ;
$dataref[ 0 ] = \$stn ;
if( defined( $id )) {
$id = padstr( $id, $staNameLen ) ;
$dataref[ 1 ] = \$id ;
$dataref[ 2 ] = \$lat ;
$dataref[ 3 ] = \$lon ;
$dataref[ 4 ] = \$elev ;
}
$dataref[ 5 ] = \theTime( "synoptic" ) ;
$dataref[ 11 ] = \theTime( "release" ) ; # release time
}
setvars() ;
$status = NetCDF::recput( $ncid, $record, \@dataref ) ;
if( $status ) { # failure
print "Wmo =$stn, Type =$type, Rpt: $rday$rhour, ",
"Trans: $bday$bhour, recput =$status\n" ;
} else {
$status = NetCDF::sync( $ncid ) ;
#print "Syncing $ncfile with ncid $ncid\n" ;
}
atexit( "eof" ) if( eof( STDIN )) ;
if( $type eq "TTAA" || $type eq "TTCC" ) {
undef( $numMand ) ;
undef( @prMan ) ;
undef( @htMan ) ;
undef( @tpMan ) ;
undef( @tdMan ) ;
undef( @wdMan ) ;
undef( @wsMan ) ;
undef( $numTrop ) ;
undef( @prTrop ) ;
undef( @tpTrop ) ;
undef( @tdTrop ) ;
undef( @wdTrop ) ;
undef( @wsTrop ) ;
undef( $numMwnd ) ;
undef( @prMaxW ) ;
undef( @wdMaxW ) ;
undef( @wsMaxW ) ;
} elsif( $type eq "TTBB" || $type eq "TTDD" ) {
undef( $numSigT ) ;
undef( @prSigT ) ;
undef( @tpSigT ) ;
undef( @tdSigT ) ;
undef( %prSigT ) ;
} elsif( $type eq "PPBB" || $type eq "PPDD" ) {
undef( $numSigW ) ;
undef( @htSigW ) ;
undef( @wsSigW ) ;
undef( @wdSigW ) ;
undef( %htSigW ) ;
}
undef( $i ) ;
undef( $j ) ;
undef( $k ) ;
undef( @pres ) ;
undef( @ht ) ;
undef( @t ) ;
undef( @td ) ;
undef( @wdir ) ;
undef( @wspd ) ;
undef( $trop_pres ) ;
undef( $trop_t ) ;
undef( $trop_td_dep ) ;
undef( $trop_wdir ) ;
undef( $trop_wspd ) ;
undef( $max_pres ) ;
undef( $max_wdir ) ;
undef( $max_wspd ) ;
} # end foreach report
} # end while( 1 )
atexit( "eof" ) ;
exit( 0 ) ; #should never get here
# setvars into record
sub setvars
{
# enter report data into proper record location
SWITCH1: {
( $type eq "TTAA" || $type eq "TTCC" ) && do {
$numMand = ${$dataref[ 6 ]} ; # number of mandatory
# data already inserted into record, append/overwrite
if( $numMand != $F ) {
@prMan = @{ $dataref[ 13 ] } ;
@htMan = @{ $dataref[ 14 ] } ;
@tpMan = @{ $dataref[ 15 ] } ;
@tdMan = @{ $dataref[ 16 ] } ;
@wdMan = @{ $dataref[ 17 ] } ;
@wsMan = @{ $dataref[ 18 ] } ;
} else { # first data write
@prMan = @ML ;
@htMan = @ML ;
@tpMan = @ML ;
@tdMan = @ML ;
@wdMan = @ML ;
@wsMan = @ML ;
}
if( $pres[ 1 ] == 1000 ) {
$prMan[ 0 ] = $pres[ 0 ] ;
$htMan[ 0 ] = $ht[ 0 ] ;
$tpMan[ 0 ] = $t[ 0 ] ;
$tdMan[ 0 ] = $td[ 0 ] ;
$wdMan[ 0 ] = $wdir[ 0 ] ;
$wsMan[ 0 ] = $wspd[ 0 ] ;
$j = 1 ;
} else {
$j = 0 ;
}
for( $i = $j; $i < $manLevel; $i++ ) {
$k = $mand_pres{ "$pres[ $i ]" } ;
next unless( defined( $k )) ;
$prMan[ $k ] = $pres[ $i ] ;
$htMan[ $k ] = $ht[ $i ] ;
$tpMan[ $k ] = $t[ $i ] ;
$tdMan[ $k ] = $td[ $i ] ;
$wdMan[ $k ] = $wdir[ $i ] ;
$wsMan[ $k ] = $wspd[ $i ] ;
}
$numMand = 0 ;
for( $i = 0 ; $i < $manLevel; $i++ ) {
$numMand++ if( $prMan[ $i ] != $F ) ;
}
$dataref[ 6 ] = \$numMand ;
$dataref[ 13 ] = \@prMan ;
$dataref[ 14 ] = \@htMan ;
$dataref[ 15 ] = \@tpMan ;
$dataref[ 16 ] = \@tdMan ;
$dataref[ 17 ] = \@wdMan ;
$dataref[ 18 ] = \@wsMan ;
# enter topopause data
$numTrop = ${$dataref[ 10 ]} ; # number of tropopause
# data already inserted into record, append/overwrite
if( $numTrop != $F ) {
@prTrop = @{ $dataref[ 25 ] } ;
@tpTrop = @{ $dataref[ 26 ] } ;
@tdTrop = @{ $dataref[ 27 ] } ;
@wdTrop = @{ $dataref[ 28 ] } ;
@wsTrop = @{ $dataref[ 29 ] } ;
} else { # first data write
@prTrop = @TL ;
@tpTrop = @TL ;
@tdTrop = @TL ;
@wdTrop = @TL ;
@wsTrop = @TL ;
}
# enter new tropopause data
if( defined( $trop_pres ) && $numTrop == $F ) {
$prTrop[ 0 ] = $trop_pres ;
$tpTrop[ 0 ] = $trop_t ;
$tdTrop[ 0 ] = $trop_td_dep ;
$wdTrop[ 0 ] = $trop_wdir ;
$wsTrop[ 0 ] = $trop_wspd ;
$numTrop = 1;
} elsif( defined( $trop_pres ) ) {
for( $i = 0; $i < $mTropNum; $i++ ) {
if( $prTrop[ $i ] == $trop_pres ) {
$index = $i ;
last ;
} elsif( $prTrop[ $i ] == $F ) {
$numTrop++ ;
$index = $i ;
last ;
}
}
$prTrop[ $index ] = $trop_pres ;
$tpTrop[ $index ] = $trop_t ;
$tdTrop[ $index ] = $trop_td_dep ;
$wdTrop[ $index ] = $trop_wdir ;
$wsTrop[ $index ] = $trop_wspd ;
}
$dataref[ 25 ] = \@prTrop ;
$dataref[ 26 ] = \@tpTrop ;
$dataref[ 27 ] = \@tdTrop ;
$dataref[ 28 ] = \@wdTrop ;
$dataref[ 29 ] = \@wsTrop ;
$dataref[ 10 ] = \$numTrop ;
# enter maximum wind data
$numMwnd = ${$dataref[ 9 ]} ; # number of max winds
# data already inserted into record, append/overwrite
if( $numMwnd != $F ) {
@prMaxW = @{ $dataref[ 30 ] } ;
@wdMaxW = @{ $dataref[ 31 ] } ;
@wsMaxW = @{ $dataref[ 32 ] } ;
} else { # first data write
@prMaxW = @MWL ;
@wdMaxW = @MWL ;
@wsMaxW = @MWL ;
}
# enter new max wind data
if( defined( $max_pres ) && $numMwnd == $F ) {
$prMaxW[ 0 ] = $max_pres ;
$wdMaxW[ 0 ] = $max_wdir ;
$wsMaxW[ 0 ] = $max_wspd ;
$numMwnd = 1;
} elsif( defined( $max_pres ) ) {
for( $i = 0; $i < $mWndNum; $i++ ) {
if( $prMaxW[ $i ] == $max_pres ) {
$index = $i ;
last ;
} elsif( $prMaxW[ $i ] == $F ) {
$numMwnd++ ;
$index = $i ;
last ;
}
}
$prMaxW[ $index ] = $max_pres ;
$wdMaxW[ $index ] = $max_wdir ;
$wsMaxW[ $index ] = $max_wspd ;
}
$dataref[ 30 ] = \@prMaxW ;
$dataref[ 31 ] = \@wdMaxW ;
$dataref[ 32 ] = \@wsMaxW ;
$dataref[ 9 ] = \$numMwnd ;
last SWITCH1 ;
} ;
( $type eq "TTBB" || $type eq "TTDD" ) && do {
$k = ${$dataref[ 7 ]} ; # number of Sig Temps
# data already inserted into record, overwrite or
# append
if( $k != $F ) {
@prSigT = @{ $dataref[ 19 ] } ;
@tpSigT = @{ $dataref[ 20 ] } ;
@tdSigT = @{ $dataref[ 21 ] } ;
for( $i = 0; $i < $k; $i++ ) {
$prSigT[ $i ] = sprintf( "%3.1f",
$prSigT[ $i ] ) ;
$prSigT{ "$prSigT[ $i ]" } = $i ;
}
for( $i = 0; $i <= $#pres; $i++ ) {
if( defined( $prSigT{ "$pres[ $i ]"} )) {
$j = $prSigT{ "$pres[ $i ]" } ;
$tpSigT[ $j ] = $t[ $i ] ;
$tdSigT[ $j ] = $td[ $i ] ;
} else { # append
$prSigT[ $k ] = $pres[ $i ] ;
$tpSigT[ $k ] = $t[ $i ] ;
$tdSigT[ $k ] = $td[ $i ] ;
$k++ ;
}
}
while( $k > $sigTLevel ) {
$k-- ;
pop( @prSigT ) ;
pop( @tpSigT ) ;
pop( @tdSigT ) ;
print "High sigTLevel $stn $type $k\n" ;
}
$dataref[ 7 ] = \$k ;
$dataref[ 19 ] = \@prSigT ;
$dataref[ 20 ] = \@tpSigT ;
$dataref[ 21 ] = \@tdSigT ;
} else { # first data entry
$numSigT = $#pres +1 ;
$dataref[ 7 ] = \$numSigT ;
for( $i = $numSigT ; $i < $sigTLevel; $i++ ) {
$pres[ $i ] = $F ;
$t[ $i ] = $F ;
$td[ $i ] = $F ;
}
$dataref[ 19 ] = \@pres ;
$dataref[ 20 ] = \@t ;
$dataref[ 21 ] = \@td ;
}
$dataref[ 12 ] = \$indicator if( $type eq "TTBB" ) ;
last SWITCH1 ;
} ;
( $type eq "PPBB" || $type eq "PPDD" ) && do {
$k = ${$dataref[ 8 ]} ; # number of Sig Winds
# data already inserted into record, overwrite or
# append
if( $k != $F ) {
@htSigW = @{ $dataref[ 22 ] } ;
@wdSigW = @{ $dataref[ 23 ] } ;
@wsSigW = @{ $dataref[ 24 ] } ;
for( $i = 0; $i < $k; $i++ ) {
$htSigW[ $i ] = sprintf( "%3.1f",
$htSigW[ $i ] ) ;
$htSigW{ "$htSigW[ $i ]" } = $i ;
}
for( $i = 0; $i <= $#ht; $i++ ) {
if( defined( $htSigW{ "$ht[ $i ]" } )) {
$j = $htSigW{ "$ht[ $i ]" } ;
$wsSigW[ $j ] = $wspd[ $i ] ;
$wdSigW[ $j ] = $wdir[ $i ] ;
} else { # append
$htSigW[ $k ] = $ht[ $i ] ;
$wsSigW[ $k ] = $wspd[ $i ] ;
$wdSigW[ $k ] = $wdir[ $i ] ;
$k++ ;
}
}
while( $k > $sigWLevel ) {
$k-- ;
pop( @htSigW ) ;
pop( @wsSigW ) ;
pop( @wdSigW ) ;
print "High sigWLevel $stn $type $k\n" ;
}
$dataref[ 8 ] = \$k ;
$dataref[ 22 ] = \@htSigW ;
$dataref[ 23 ] = \@wdSigW ;
$dataref[ 24 ] = \@wsSigW ;
} else { # first data entry
$numSigW = $#wdir +1 ;
$dataref[ 8 ] = \$numSigW ;
for( $i = $numSigW; $i < $sigWLevel; $i++ ) {
$ht[ $i ] = $F ;
$wdir[ $i ] = $F ;
$wspd[ $i ] = $F ;
}
$dataref[ 22 ] = \@ht ;
$dataref[ 23 ] = \@wdir ;
$dataref[ 24 ] = \@wspd ;
}
last SWITCH1 ;
} ;
}
} # end setvars
# create a netcdf file or reopen a existing one
sub doNet
{
my( $Ncfile, $Id, $Num, $Time, $baseTime, $offset, $rpt, $stn ) ;
#$ncfile = $datadir . "/" . $yyyymmddhh . "_perl_ua.nc" ;
$ncfile = $datadir . "/" . $upperairPrefix . substr( $yyyymmddhh, 0, 8 )
. $upperair . ".nc" ;
if( $hourly ) {
$ncfile = $datadir . "/" . $upperairPrefix . substr( $yyyymmddhh, 0, 8
) .
"_" . substr( $yyyymmddhh, 8 ) . $upperair . ".nc" ;
}
# writing to same file
return 1 if( $ncfile eq $lastNc ) ;
# current time
$thetime = time() ;
# save current file info
$Nets{ $lastNc } = "$ncid $recnum $thetime" if( $lastNc ) ;
# File is open, get ncfile id and recnum and reset the time
if( defined( $Nets{ $ncfile } )) { # already open for writes
( $ncid, $recnum, $ncTime ) = split( " ", $Nets{ $ncfile } ) ;
$ncTime = $thetime ;
$lastNc = $ncfile ;
return 1 ;
}
# close files with no activity for 20 minutes
foreach $Ncfile ( keys %Nets ) {
( $Id, $Num, $Time ) = split( " ", $Nets{ $Ncfile } ) ;
if( $thetime - $Time > 1200 ) {
print "Closing $Ncfile with ncid $Id, No write for > 20
Minutes\n" ;
$status = NetCDF::close( $Id ) ;
delete( $Nets{ $Ncfile } ) ;
# } elsif( $Ncfile eq $lastNc ) {
# print "Syncing $lastNc with ncid $ncid\n" ;
# $status = NetCDF::sync( $ncid ) ;
}
}
# remove rpt entries older than 24 hours
$baseTime = $yyyymm . $cday . $chour ;
foreach $rpt ( keys %rpt_hash ) {
( $stn, $Time ) = split( " ", $rpt ) ;
$offset = $baseTime - $Time ;
next if( $offset < 100 ) ; # same day ok & previous day ok > $chour
delete( $rpt_hash{ $rpt } ) ;
}
# open or create ncfiles
if( -e $ncfile ) {
$ncid = NetCDF::open( "$ncfile", WRITE ) ;
return 0 if( $ncid == -1 ) ;
$recNum_id = NetCDF::dimid( $ncid, "recNum" ) ;
$name_id = "xxxxxxxx" ;
$recnum = -1 ;
# get current value of recnum
NetCDF::diminq( $ncid, $recNum_id, $name_id, $recnum ) ;
} else {
system( "$ncgen -o $ncfile $cdlfile" ) ;
$ncid = NetCDF::open( "$ncfile", WRITE ) ;
return 0 if( $ncid == -1 ) ;
# NetCDF record counter
$recnum = 0 ;
}
$Nets{ $ncfile } = "$ncid $recnum $thetime" ;
$lastNc = $ncfile ;
print "Opening $ncfile with ncid $ncid\n" ;
return 1 ;
} # end doNet
# print vars
sub printvars
{
my( @rec ) ;
for( $i=0; $i<$j; $i++ ) {
$rec[$i] = sprintf "%6s%8s%8s%10s%10s%13s",
$pres[$i], $ht[$i], $t[$i], $td[$i], $wdir[$i], $wspd[$i] ;
}
# Output decoded sounding
print "Station: $stn, report type: $type, j = $j\n" ;
print "Time: $yyyymmddhh\n" ;
print "---------------------------------------------------------\n" ;
print "Pres Ht Temp DewPt Dir Spd\n" ;
print " mb m C C deg m/s \n" ;
print "---------------------------------------------------------\n" ;
for( $i = 0; $i < $j; $i++ ) { print "$rec[$i]\n" ; }
print "\n" ;
if( $trop_t ) {
print "Tropopause values\n" ;
print "---------------------------------------------------------\n" ;
# trop_t, trop_td_dep trop_wdir, trop_wspd
$trop = sprintf "%6s%15s%10s%10s%13s",
$trop_pres, $trop_t, $trop_td_dep, $trop_wdir, $trop_wspd ;
print "$trop\n\n" ;
}
if( $max_wdir ) {
# max_pres, max_wdir, max_wspd
print "Maximum values\n" ;
print "---------------------------------------------------------\n" ;
$max = sprintf "%6s%35s%13s", $max_pres, $max_wdir, $max_wspd ;
print "$max\n\n" ;
}
} # end printvars
# execute at exit
sub atexit
{
my( $sig ) = @_ ;
if( $sig eq "eof" ) {
print "eof on STDIN --shutting down\n" ;
} elsif( defined( $sig )) {
print "Caught SIG$sig --shutting down\n" ;
}
# open ua.lst, list of reports processed in the last 8 hours.
open( LST, ">$datadir/ua.lst" ) ||
die "could not open $datadir/ua.lst: $!\n" ;
select( LST ) ;
# remove rpt entries older than 24 hours
$yyyymmddhh = $yyyymm . $cday . $chour ;
foreach $rpt ( keys %rpt_hash ) {
( $stn, $Time ) = split( " ", $rpt ) ;
$offset = $yyyymmddhh - $Time ;
next unless( $offset < 100 ) ; # same day ok & previous day ok > $chour
print "$rpt $rpt_hash{ $rpt }\n" ;
}
close LST ;
foreach $Ncfile ( keys %Nets ) {
( $ncid, $recnum, $nctime ) = split( " ", $Nets{ $Ncfile } ) ;
print STDOUT "Closing $Ncfile with ncid $ncid\n" ;
$status = NetCDF::close( $ncid ) ;
print "STDOUT NetCDF::close status = $status\n" if( $status ) ;
}
close( STDOUT ) ;
close( STDERR ) ;
exit( 0 ) ;
} # end atexit
# pad str to correct length
sub padstr
{
( $str, $len ) = @_ ;
my( $size, $i ) ;
$size = length( $str ) ;
for( $i = $size; $i < $len; $i++ ) {
$str .= "\0" ;
#print "$str,\n" ;
}
if( $size > $len ) {
print STDOUT "String length is over $len chars long:\n $str\n"
if( $verbose ) ;
$str = substr( $str, 0, $len ) ;
}
return $str ;
} # end padstr
sub makedataref
{
undef( @dataref ) ;
$F = 99999 ;
$F0 = 99999 ;
$F1 = 99999 ;
$F2 = 99999 ;
$F3 = 99999 ;
$F4 = 99999 ;
$F5 = 99999 ;
$F6 = 99999 ;
$F7 = 99999 ;
$F8 = 99999 ;
$F85 = 99999 ;
$F9 = 99999 ;
$F10 = 99999 ;
$A = \$F ;
$S6 = "\0" x $staNameLen ;
$AS6 = \$S6 ;
@ML = ( 99999 ) x ( $manLevel ) ;
@ML0 = ( 99999 ) x ( $manLevel ) ;
@ML1 = ( 99999 ) x ( $manLevel ) ;
@ML2 = ( 99999 ) x ( $manLevel ) ;
@ML3 = ( 99999 ) x ( $manLevel ) ;
@ML4 = ( 99999 ) x ( $manLevel ) ;
@ML5 = ( 99999 ) x ( $manLevel ) ;
@STL = ( 99999 ) x ( $sigTLevel ) ;
@STL0 = ( 99999 ) x ( $sigTLevel ) ;
@STL1 = ( 99999 ) x ( $sigTLevel ) ;
@STL2 = ( 99999 ) x ( $sigTLevel ) ;
@SWL = ( 99999 ) x ( $sigWLevel ) ;
@SWL0 = ( 99999 ) x ( $sigWLevel ) ;
@SWL1 = ( 99999 ) x ( $sigWLevel ) ;
@SWL2 = ( 99999 ) x ( $sigWLevel ) ;
@TL = ( 99999 ) x ( $mTropNum ) ;
@TL0 = ( 99999 ) x ( $mTropNum ) ;
@TL1 = ( 99999 ) x ( $mTropNum ) ;
@TL2 = ( 99999 ) x ( $mTropNum ) ;
@TL3 = ( 99999 ) x ( $mTropNum ) ;
@TL4 = ( 99999 ) x ( $mTropNum ) ;
@MWL = ( 99999 ) x ( $mWndNum ) ;
@MWL0 = ( 99999 ) x ( $mWndNum ) ;
@MWL1 = ( 99999 ) x ( $mWndNum ) ;
@MWL2 = ( 99999 ) x ( $mWndNum ) ;
# default netCDF record structure, contains all vars for the UA reports
@dataref = ( \$F0, \$S6, \$F1, \$F2, \$F3, \$F4, \$F5, \$F6, \$F7, \$F8,
\$F85, \$F9, \$F10, \@ML0, \@ML1, \@ML2, \@ML3, \@ML4, \@ML5, \@STL0, \@STL1,
\@STL2, \@SWL0, \@SWL1, \@SWL2, \@TL0, \@TL1, \@TL2, \@TL3, \@TL4, \@MWL0,
\@MWL1, \@MWL2 ) ;
}
sub theTime
{
my( $ss, $mm, $hh, $mday, $mon, $fmon, $year, $fyear, $wday, $yday, $isdst ) ;
my( $when ) = @_ ;
$mm = 0 ;
if( $when eq "synoptic" ) {
$mday = substr( $yyyymmddhh, 6, 2 ) ;
$hh = substr( $yyyymmddhh, 8, 2 ) ;
$fmon = substr( $yyyymmddhh, 4, 2 ) ;
$fyear = substr( $yyyymmddhh, 2, 2 ) ;
} elsif( $when eq "release" ) {
$mday = $rday ;
$hh = $rhour ;
$fmon = $themonth ;
$fyear = $thedecade ;
}
$time = timegm( 0, $mm, $hh, $mday, $fmon -1, $fyear, 0, 0, 0 ) ;
return $time ;
} # end theTime
# TTAA type of reports, mandatory levels: surface to 100 mb
sub ttaa
{
my( $TTAA, $stride, $toplevel, $offset, $thermo, $wind ) ;
( $TTAA ) = @_ ;
$stride = 18 ; # when winds present in report
$toplevel = $TTAA_top_wind_level{ "$indicator" } ;
# mandatory levels, sfc thru 100 mb
for( $offset = 0 ; $offset < length( $TTAA ) ; $offset += $stride, $j++ ) {
$pres[ $j ] = $ht[ $j ] = $t[ $j ] = $td[ $j ] = $wdir[ $j ] =
$wspd[ $j ] = $F ;
$pres[ $j ] = substr( $TTAA, $offset, 2 ) ;
last unless( defined( $mand_ttaa{ "$pres[ $j ]" } ) ) ;
$pres[ $j ] = $mand_ttaa{ "$pres[ $j ]" } ;
# $stride is set to 12 when no winds are reported
$stride = 12 if( $pres[ $j ] < $toplevel ) ;
if( substr( $TTAA, $offset+2, 10 ) =~ m#([\d/]{3}) ([\d/]{5}) # ) {
$ht[ $j ] = $1 ;
dec_thermo( $2 ) ;
} else {
$ht[ $j ] = "///" ;
}
if(( $stride == 18 || $pres[ $j ] == 99999 ) &&
substr( $TTAA, $offset + 12, 6 ) =~ m#([\d/]{5}) # ) {
dec_wind( $1 ) ;
}
if( $pres[ $j ] == 99999 ) { # surface
$sfc_pres = $pres[ $j ] =
( $ht[ $j ] < 100 ) ? $ht[ $j ] + 1000 : $ht[ $j ] ;
$ht[ $j ] = $elev ;
}
$ht[ $j ] .= "0" if( $pres[ $j ] <= 500 ) ; # 500 mb and above
if( $ht[ $j ] =~ m#/# ) {
$ht[ $j ] = $F ;
# 1000 mb special case for below sea level
} elsif( $pres[ $j ] == 1000 && $pres[ $j ] > $sfc_pres ) {
if( $ht[ $j ] > 500 ) {
$ht[ $j ] -= 500 ;
$ht[ $j ] *= -1 ;
}
# 925 mb special case for below sea level
} elsif( $pres[ $j ] == 925 && $pres[ $j ] > $sfc_pres ) {
if( $ht[ $j ] > 500 ) {
$ht[ $j ] -= 500 ;
$ht[ $j ] *= -1 ;
}
} elsif( $pres[ $j ] == 925 && $ht[ $j ] < 250 ) { # 925 mb
$ht[ $j ] = "1" . $ht[ $j ] ;
} elsif( $pres[ $j ] == 850 ) { # 850 mb
$ht[ $j ] = "1" . $ht[ $j ] ;
} elsif( $pres[ $j ] == 700 ) { # 700 mb
if( $ht[ $j ] < 350 ) {
$ht[ $j ] = "3" . $ht[ $j ] ;
} else {
$ht[ $j ] = "2" . $ht[ $j ] ;
}
} elsif( $pres[ $j ] == 250 && $ht[ $j ] < 5000 ) { # 250 mb
$ht[ $j ] = "1" . $ht[ $j ] ;
} elsif( $pres[ $j ] =~ /200|150|100$/ ) { # 200, 150, &100 mb
$ht[ $j ] = "1" . $ht[ $j ] ;
}
} # End of loop thru TTAA mandatory levels
pop( @pres ) unless( $mand_pres{ "$pres[ $j ]" } ) ;
# Tropopause near end of TTAA
if( $TTAA !~ / 88999/ &&
$TTAA =~ m# 88(\d\d\d) ([\d/]{5}) ([\d/]{5}) # ) {
$trop_pres = $1 ;
$thermo = $2 ;
$wind = $3 ;
# sets trop_t, trop_td_dep
trop_dec_thermo( $thermo ) ;
# sets trop_wdir, trop_wspd
trop_dec_wind( $wind ) ;
} # End of tropopause level in TTAA group
# Maximum Wind
if( $TTAA !~ / 77999/ &&
$TTAA =~ m# (66|77)(\d\d\d) ([\d/]{5}) # ) {
$max_pres = $2 ;
# sets max_wdir, max_wspd
max_dec_wind( $3 ) ;
} # End of maximum wind level in TTAA group
} #end ttaa
# TTBB reports, significant levels with respect to temp: surface to 100 mb
sub ttbb
{
my( $TTBB, $offset ) ;
( $TTBB ) = @_ ;
$indicator = $F if( $indicator eq '/' ) ;
$offset = 0 ;
while( $offset < length( $TTBB )) {
last unless substr( $TTBB, $offset, 2 ) =~ /00|11|22|33|44|55|66|77|88/;
$pres[ $j ] = $ht[ $j ] = $t[ $j ] = $td[ $j ] = $wdir[ $j ] =
$wspd[ $j ] = $F ;
if( substr( $TTBB, $offset+2, 10 ) =~ /(\d{3}) (\d{5}) / ) {
$pres[ $j ] = $1 ;
dec_thermo( $2 ) ;
$pres[ $j ] = "1" . $pres[ $j ] if( $pres[ $j ] < 100 ) ;
if( defined( $mand_pres{ "$pres[ $j ]" } )) {
$offset += 12 ;
next ;
}
$pres[ $j ] = sprintf( "%3.1f", $pres[ $j ] ) ;
$j++ ;
}
$offset += 12 ;
} # End of TTBB
} #end ttbb
# TTCC type of reports, mandatory levels: under 100 mb
sub ttcc
{
my( $TTCC, $stride, $toplevel, $offset, $thermo, $wind ) ;
( $TTCC ) = @_ ;
$stride = 18 ; # when winds present in report
$toplevel = $TTCC_top_wind_level{ "$indicator" } ;
if( ! defined( $toplevel )) {
print "Undefined Wind termination level: $level, $stn $type ",
"Rpt: $rday$rhour Trans: $bday$bhour\n" ;
return ;
}
# mandatory levels, 70 mb thru 1 mb
for( $offset = 0 ; $offset < length( $TTCC ) ; $offset += $stride, $j++ ) {
$pres[ $j ] = $ht[ $j ] = $t[ $j ] = $td[ $j ] = $wdir[ $j ] =
$wspd[ $j ] = $F ;
$pres[ $j ] = substr( $TTCC, $offset, 2 ) ;
last unless( defined( $mand_ttcc{ "$pres[ $j ]" } ) ) ;
# $stride is set to 12 when no winds are reported
$stride = 12 if( $pres[ $j ] < $toplevel ) ;
if( substr( $TTCC, $offset+2, 10 ) =~ m#([\d/]{3}) ([\d/]{5}) # ) {
$ht[ $j ] = $1 ;
dec_thermo( $2 ) ;
} else {
$ht[ $j ] = "///" ;
}
if( $stride == 18 &&
substr( $TTCC, $offset +12, 6 ) =~ m#([0-3/][\d/]{4}) # ) {
dec_wind( $1 ) ;
}
if( $ht[ $j ] =~ m#/# ) {
$ht[ $j ] = $F ;
} else {
$ht[ $j ] .= "0" ;
}
if( $pres[ $j ] == 70 ) { # 70 mb
$ht[ $j ] = "1" . $ht[ $j ] ;
} elsif( $pres[ $j ] =~ /50|30|20/ ) { # 50,30,20 mb
$ht[ $j ] = "2" . $ht[ $j ] ;
} else { # 10, 7, 5, 3, 2, 1 mb
$ht[ $j ] = "3" . $ht[ $j ] ;
}
} # End of loop thru TTCC mandatory levels
pop( @pres ) unless( $mand_pres{ "$pres[ $j ]" } ) ;
# Tropopause near end of TTCC if not found in TTAA group
if( $TTCC !~ / 88999/ && $TTCC =~ m# 88(\d\d)(\d) ([\d/]{5}) ([\d/]{5}) # ) {
$trop_pres = $1 . "." . $2 ;
$thermo = $3 ;
$wind = $4 ;
$trop_pres += 0.5 ;
$trop_pres =~ s#\.\d*## ;
# sets trop_t, trop_td_dep
trop_dec_thermo( $thermo ) ;
# sets trop_wdir, trop_wspd
trop_dec_wind( $wind ) ;
} # End of tropopause level in TTCC group
# Maximum Wind
if( $TTCC !~ / 77999/ && $TTCC =~ m# (66|77)(\d\d)(\d) ([\d/]{5}) # ) {
$max_pres = $2 . "." . $3 ;
max_dec_wind( $4 ) ;
$max_pres += 0.5 ;
$max_pres =~ s#\.\d*## ;
}# End of maximum wind level in TTCC group
} #end ttcc
# TTDD reports, significant levels with respect to temperature: under 100 mb
sub ttdd
{
my( $TTDD, $offset ) ;
( $TTDD ) = @_ ;
$offset = 0 ;
while( $offset < length( $TTDD )) {
last unless substr( $TTDD, $offset, 2 ) =~ /00|11|22|33|44|55|66|77|88/;
$pres[ $j ] = $ht[ $j ] = $t[ $j ] = $td[ $j ] = $wdir[ $j ] =
$wspd[ $j ] = $F ;
if( substr( $TTDD, $offset +2, 10 ) =~ /(\d{3}) (\d{5}) / ) {
$pres[ $j ] = sprintf( "%3.1f", $1 / 10 ) ;
dec_thermo( $2 ) ;
$j++ ;
}
$offset += 12 ;
}
} #end ttdd
# PP(BB|DD) reports, significant levels with respect to winds: all levels
sub ppbbdd
{
my( $PP, $offset, $i, $ten_thousand, $thousand ) ;
( $PP ) = @_ ;
$offset = 0 ;
while( $offset < length( $PP )) {
last unless substr( $PP, $offset, 1 ) eq "9" ;
$ten_thousand = substr( $PP, $offset +1, 1 ) ;
for( $i = 0; $i < 3; $i++ ) {
last if(( $thousand =
substr( $PP, $offset +2 + $i, 1 )) eq "/" ) ;
$pres[ $j ] = $ht[ $j ] = $t[ $j ] = $td[ $j ] = $wdir[ $j ] =
$wspd[ $j ] = $F ;
$ht[ $j ] = ( $ten_thousand . $thousand . "000" ) * 0.3048 ;
next unless $ht[ $j ] ;
if( substr( $PP, $offset +6 + $i * 6, 6 ) =~ /(\d{5}) / ) {
dec_wind( $1 ) ;
}
next if( $wdir[ $j ] eq $F || $wspd[ $j ] eq $F ) ;
$ht[ $j ] = sprintf( "%3.1f", $ht[ $j ] ) ;
$j++ ;
}
$offset += 6 + $i * 6 ;
}
} # end ppbbdd
sub dec_thermo
{
my( $thermo ) = @_ ;
my( $td_dep ) ;
$t[ $j ] = substr( $thermo, 0, 2 ) . "." . substr( $thermo, 2, 1 ) ;
$td_dep = substr( $thermo, 3, 2 ) ;
if( $t[ $j ] =~ m#/# ) {
$t[ $j ] = $F ;
} else {
$t[ $j ] =~ s/^0// ;
if( $t[ $j ] =~ m#[13579]$# ) { $t[ $j ] = "-" . $t[ $j ] ; }
}
if( $td_dep =~ m#/# ) {
$td[ $j ] = $F ;
} else {
# change to report Dew-point depression instead of Dew-point
if( $td_dep <= 50 ) {
#$td[ $j ] = $t[ $j ] - $td_dep / 10.0 ;
$td[ $j ] = $td_dep / 10.0 ;
} else {
#$td[ $j ] = $t[ $j ] - ( $td_dep -50.0 ) ;
$td[ $j ] = $td_dep -50.0 ;
}
$td[ $j ] *= 1.0 ;
}
# Convert to kelvin
$t[ $j ] = celsius2kelvin( $t[ $j ] ) ;
return ;
} # End dec_thermo
sub trop_dec_thermo
{
my( $thermo ) = @_ ;
my( $td_dep ) ;
$trop_t = substr( $thermo, 0, 2 ) . "." . substr( $thermo, 2, 1 ) ;
$td_dep = substr( $thermo, 3, 2 ) ;
if( $trop_t =~ m#/# ) {
$trop_t = $F ;
} else {
$trop_t =~ s/^0// ;
$trop_t = "-" . $trop_t if( $trop_t =~ m#[13579]$# ) ;
}
if( $td_dep =~ m#/# ) {
$trop_td_dep = $F ;
} else {
# change to report Dew-point depression instead of Dew-point
if( $td_dep <= 50 ) {
#$trop_td_dep = $trop_t - $td_dep/10.0 ;
$trop_td_dep = $td_dep/10.0 ;
} else {
#$trop_td_dep = $trop_t - ( $td_dep -50.0 ) ;
$trop_td_dep = $td_dep -50.0 ;
}
$trop_td_dep *= 1.0 ;
}
# Convert to kelvin
$trop_t = celsius2kelvin( $trop_t ) ;
} # end trop_dec_thermo
#
sub dec_wind
{
my( $wind ) = @_ ;
$wdir[ $j ] = substr( $wind, 0, 2 ) ;
$wspd[ $j ] = substr( $wind, 2, 3 ) ;
if( $wdir[ $j ] =~ m#/# ) {
$wdir[ $j ] = $wspd[ $j ] = $F ;
} else {
$wdir[ $j ] =~ s/^0// ;
$wdir[ $j ] .= "0" ;
if( $wspd[ $j ] > 500 ) {
$wdir[ $j ] += 5 ;
$wspd[ $j ] -= 500 ;
}
}
} # end dec_wind
#
sub trop_dec_wind
{
my( $wind ) = @_ ;
$trop_wdir = substr( $wind, 0, 2 ) ;
$trop_wspd = substr( $wind, 2, 3 ) ;
if( $trop_wdir =~ m#/# ) {
$trop_wdir = $trop_wspd = $F ;
} else {
$trop_wdir =~ s/^0// ;
$trop_wdir .= "0" ;
if( $trop_wspd > 500 ) {
$trop_wdir += 5 ;
$trop_wspd -= 500 ;
}
$trop_wspd *= 0.5144 if( $knots ) ;
}
} # end trop_dec_wind
# calculate max winds
sub max_dec_wind
{
my( $wind ) = @_ ;
$max_wdir = substr( $wind, 0, 2 ) ;
$max_wspd = substr( $wind, 2, 3 ) ;
if( $max_wdir =~ m#/# ) {
$max_wdir = $max_wspd = $F ;
} else {
$max_wdir =~ s/^0// ;
$max_wdir .= "0" ;
if( $max_wspd > 500 ) {
$max_wdir += 5 ;
$max_wspd -= 500 ;
}
$max_wspd *= 0.5144 if( $knots ) ;
}
} # end max_dec_wind
# returns kelvin given celsius temperature
sub celsius2kelvin
{
my( $celsius ) = @_ ;
if( $celsius == $F ) {
return $F ;
} else {
return $celsius + 273.15 ;
}
} # end celsius2kelvin