c Subroutine TRUEWIND c Converted from IDL true wind computation code true_wind.pro. c Created: 12/17/96 c Comments updated: 10/01/97 c Developed by: Shawn R. Smith and Mark A. Bourassa c Programmed by: Mylene Remigio c Last updated: 9/30/2014 c Direct questions to: wocemet@coaps.fsu.edu c 9/30/2014 :Version 2 - Arturo Valery c If the true wind has speed and its coming from the north c then its direction should be 360deg. The problem was fixed c in which the programs’s output showed 0deg instead of 360deg. c This routine will compute meteorological true winds (direction from c which wind is blowing, relative to true north; and speed relative to c the fixed earth). c INPUT VALUES: c num integer Number of observations in input (crse, cspd, c wdir, wspd, hd) and output (adir, tdir, tspd) c data arrays. ALL ARRAYS MUST BE OF EQUAL LENGTH. c sel integer Sets option for diagnostic output. There are c four settings: c c Option 4: Calculates true winds from input c arrays with no diagnostic output or c warnings. NOT RECOMMENDED. c Option 3: [DEFAULT] Diagnostic output lists the c array index and corresponding c variables that either violate the c range checks or are equal to the c missing value. An additional table c lists the number of observation times c with no missing values, some (but not c all) missing values, and all missing c values; as well as similar totals for c the observation times that fail the c range checks. Range checks identify c negative input values and verify c directions to be between 0 and 360 c degrees. c Option 2: In addition to the default c diagnostics (option 3), a table of c all input and output values for c observation times with missing data c is provided. c Option 1: Full diagnostics -- In addition to c the diagnostics provided by option 2 c and 3, a complete data chart is c output. The table contains input and c output values for all observation c times passed to truewind. c crse real array Course TOWARD WHICH the vessel is moving over c the ground. Referenced to true north and the c fixed earth. c cspd real array Speed of vessel over the ground. Referenced c to the fixed earth. c hd real array Heading toward which bow of vessel is pointing. c Referenced to true north. c zlr real Zero line reference -- angle between bow and c zero line on anemometer. Direction is clockwise c from the bow. (Use bow=0 degrees as default c when reference not known.) c wdir real array Wind direction measured by anemometer, c referenced to the ship. c wspd real array Wind speed measured by anemometer,referenced to c the vessel's frame of reference. c wmis real array Five element array containing missing values for c crse, cspd, wdir, wspd, and hd. In the output, c the missing value for tdir is identical to the c missing value specified in wmis for wdir. c Similarly, tspd uses the missing value assigned c to wmis for wspd. c *** WDIR MUST BE METEOROLOGICAL (DIRECTION FROM)! CRSE AND CSPD MUST c BE RELATIVE TO A FIXED EARTH! *** c OUTPUT VALUES: c tdir real array True wind direction - referenced to true north c and the fixed earth with a direction from which c the wind is blowing (meteorological). c tspd real array True wind speed - referenced to the fixed earth. c adir real array Apparent wind direction (direction measured by c wind vane, relative to true north). IS c REFERENCED TO TRUE NORTH & IS DIRECTION FROM c WHICH THE WIND IS BLOWING. Apparent wind c direction is the sum of the ship relative wind c direction (measured by wind vane relative to the c bow), the ship's heading, and the zero-line c reference angle. NOTE: The apparent wind speed c has a magnitude equal to the wind speed measured c by the anemometer. c DIAGNOSTIC OUTPUT: c nw integer Number of observation times for which tdir and c tspd were calculated (without missing values) c nwpm integer Number of observation times with some values c (crse, cspd, wdir, wspd, hd) missing. tdir, c tspd set to missing value. c nwam integer Number of observation times with all values c (crse, cspd, wdir, wspd, hd) missing. tdir, c tspd set to missing value. c nwf integer Number of observation times where the program c fails -- at least one of the values (crse, cspd, c wdir, wspd, hd) is invalid ************************************************************************ subroutine truewind (num,sel,crse,cspd,wdir,zlr,hd,adir,wspd,wmis, + tdir,tspd,nw,nwpm,nwam,nwf) c DEFINE VARIABLES real crse(num), cspd(num), hd(num), wdir(num), adir(num), zlr, + wspd(num), wmis(5), tdir(num), tspd(num), x, y, mcrse, mwdir, + mtdir, dtor, pi integer num, nw, nwpm, nwam, nwf, calm_flag, sel c INITIALIZE VALUES. x = 0.0 y = 0.0 nw = 0 nwpm = 0 nwam = 0 nwf = 0 pi = ACOS(-1.) dtor = pi/180. ! DEGREES TO RADIAN CONVERSION c LOOP OVER 'NUM' VALUES. do i = 1,num c CHECK COURSE, SHIP SPEED, HEADING, WIND DIRECTION, AND WIND c SPEED FOR VALID VALUES (i.e. neither missing nor outside c physically acceptable ranges). 30 if( ( ((crse(i) .lt. 0.0).or.(crse(i) .gt. 360.0)) .and. + (crse(i) .ne. wmis(1)) ) .or. ( (cspd(i) .lt. 0.0).and. + (cspd(i) .ne. wmis(2)) ) .or. ( ((wdir(i) .lt. 0.0).or. + (wdir(i) .gt. 360.0)) .and. (wdir(i) .ne. wmis(3)) ) .or. + ( (wspd(i) .lt. 0.0).and.(wspd(i) .ne. wmis(4)) ) .or. + ( ( (hd(i) .lt. 0.0).or.(hd(i) .gt. 360.0)) .and. + (hd(i) .ne. wmis(5)) ) )then C WHEN SOME OR ALL OF INPUT DATA FAILS RANGE CHECK, SET TRUE C WINDS TO MISSING. STEP INDEX FOR INPUT VALUE(S) BEING OUT C OF RANGE nwf = nwf + 1 tdir(i) = wmis(3) tspd(i) = wmis(4) if( (crse(i) .ne. wmis(1)) .and. (cspd(i) .ne. wmis(2)) .and. + (wdir(i) .ne. wmis(3)) .and. (wspd(i) .ne. wmis(4)) .and. + (hd(i) .ne. wmis(5)) )then C STEP INDEX FOR ALL INPUT VALUES BEING NON-MISSING nw = nw + 1 elseif( (crse(i) .ne. wmis(1)) .or. (cspd(i) .ne. wmis(2)) + .or. (wdir(i) .ne. wmis(3)) .or. (wspd(i) .ne. wmis(4)) + .or. (hd(i) .ne. wmis(5)) )then C STEP INDEX FOR PART OF INPUT VALUES BEING MISSING nwpm = nwpm + 1 else C STEP INDEX FOR ALL INPUT VALUES BEING MISSING nwam = nwam + 1 endif c WHEN COURSE, SHIP SPEED, HEADING, WIND DIRECTION, AND WIND c SPEED ARE ALL IN RANGE AND NON-MISSING, THEN COMPUTE TRUE WINDS. elseif( (crse(i) .ne. wmis(1)) .and. (cspd(i) .ne. wmis(2)) + .and. (wdir(i) .ne. wmis(3)) .and. (wspd(i) .ne. wmis(4)) + .and. (hd(i) .ne. wmis(5)) )then nw = nw +1 C CONVERT FROM NAVIGATIONAL COORDINATES TO ANGLES COMMONLY USED C IN MATHEMATICS mcrse = 90.0 - crse(i) if(mcrse .le. 0.0) + mcrse = mcrse + 360.0 ! KEEP VALUES BETWEEN 0 AND 360 DEGREES C CHECK ZLR FOR VALID VALUE. IF NOT VALID, SET EQUAL TO ZERO. if((zlr .lt. 0.0) .or. (zlr .gt. 360)) + zlr = 0.0 C CALCULATE APPARENT WIND DIRECTION adir(i) = hd(i) + wdir(i) + zlr C KEEP ADIR BETWEEN 0 AND 360 DEGREES do while (adir(i) .ge. 360.0) adir(i) = adir(i) - 360 enddo C CONVERT FROM METEOROLOGICAL COORDINATES TO ANGLES COMMONLY C USED IN MATHEMATICS mwdir = 270.0 - adir(i) C KEEP MWDIR BETWEEN 0 AND 360 DEGREES if (mwdir .le. 0.0) + mwdir = mwdir + 360.0 if (mwdir .gt. 360.0) + mwdir = mwdir - 360.0 C DETERMINE THE EAST-WEST AND NORTH-SOUTH VECTOR COMPONENTS OF C THE TRUE WIND x = wspd(i) * COS(mwdir * dtor) + cspd(i) * COS(mcrse * dtor) y = wspd(i) * SIN(mwdir * dtor) + cspd(i) * SIN(mcrse * dtor) C USE THE TWO VECTOR COMPONENTS TO CALCULATE THE TRUE WIND SPEED tspd(i) = SQRT(x * x + y * y) calm_flag = 1 C DETERMINE THE ANGLE FOR THE TRUE WIND if(ABS(x) .gt. 0.00001) then mtdir = ATAN2(y,x) / dtor else if(ABS(y) .gt. 0.00001) then mtdir = 180.0 - 90.0 * (y / ABS(y)) else mtdir = 270.0 calm_flag = 0 endif endif C CONVERT FROM THE COMMON MATHEMATICAL ANGLE COORDINATE TO THE C METEOROLOGICAL WIND DIRECTION tdir(i) = 270.0 - mtdir C MAKE SURE THAT THE TRUE WIND DIRECTION IS BETWEEN 0 AND 360 C DEGREES 10 if(tdir(i) .lt. 0.0) then tdir(i) = (tdir(i) + 360.0) * calm_flag goto 10 endif 20 if(tdir(i) .gt. 360.0) then tdir(i) = (tdir(i) - 360.0) * calm_flag goto 20 endif if ((calm_flag .eq. 1) .and. (tdir(i) .lt. 0.0001)) then tdir(i) = 360.0 endif x = 0.0 y = 0.0 c WHEN COURSE, SHIP SPEED, APPARENT DIRECTION, AND WIND SPEED c ARE ALL IN RANGE BUT PART OF THESE INPUT VALUES ARE MISSING, C THEN SET TRUE WIND DIRECTION AND SPEED TO MISSING. elseif( (crse(i) .ne. wmis(1)) .or. (cspd(i) .ne. wmis(2)) .or. + (wdir(i) .ne. wmis(3)) .or. (wspd(i) .ne. wmis(4)) .or. + (hd(i) .ne. wmis(5)) ) then nwpm = nwpm + 1 tdir(i) = wmis(3) tspd(i) = wmis(4) c WHEN COURSE, SHIP SPEED, APPARENT DIRECTION, AND WIND SPEED c ARE ALL IN RANGE BUT ALL OF THESE INPUT VALUES ARE MISSING, C THEN SET TRUE WIND DIRECTION AND SPEED TO MISSING. else nwam = nwam + 1 tdir(i) = wmis(3) tspd(i) = wmis(4) endif enddo c OUTPUT OPTION SELECTION FOR TRUEWIND.F if(sel .eq. 1) then call full(num,crse,cspd,wdir,zlr,hd,adir,wspd,tdir,tspd) call missing_values(num,crse,cspd,wdir,hd,wspd,tdir,tspd,wmis) call truerr(num,crse,cspd,hd,wdir,wspd,wmis,nw,nwpm,nwam,nwf) else if(sel .eq. 2) then call missing_values(num,crse,cspd,wdir,hd,wspd,tdir,tspd,wmis) call truerr(num,crse,cspd,hd,wdir,wspd,wmis,nw,nwpm,nwam,nwf) else if(sel .eq. 3) then call truerr(num,crse,cspd,hd,wdir,wspd,wmis,nw,nwpm,nwam,nwf) else if(sel .eq. 4) then continue else print*,'Number invalid; using sel = 3 as default.' call truerr(num,crse,cspd,hd,wdir,wspd,wmis,nw,nwpm,nwam,nwf) endif return end C ********************************************************************** c OUTPUT SUBROUTINES C ********************************************************************** c Subroutine: FULL.F c Purpose: Produces a complete data table with all values. c Accessed only when selection #1 is chosen. subroutine full(num,crse,cspd,wdir,zlr,hd,adir,wspd,tdir,tspd) real crse(num), cspd(num), wdir(num), zlr, hd(num), adir(num), + wspd(num), tdir(num), tspd(num) integer num character cv*22, wd*24, vs*22, ws*22, str1*30, str2*10, str3*25, + str4*10 print*,'---------------------------------------------------------- +-------------------------' print*,'\n FULL TABLE' print*,' ************' print*,' index course sspeed windir zeroln shiphd | appspd | + appdir trudir truspd' do i=1,num write(*,50) i, crse(i), cspd(i), wdir(i), zlr, hd(i), wspd(i), + adir(i), tdir(i), tspd(i) 50 format(i8, 5f8.1, ' | ', f6.1, ' |', 3f8.1) enddo print*,' ' print*,' Note: Wind speed measured by anemometer is +identical' print*,' to apparent wind speed (appspd).' return end C ********************************************************************** c Subroutine: MISSING_VALUES.F c Purpose: Produces a data table of the data with missing values. c Accessed when selection #1 or #2 is chosen. subroutine missing_values(num,crse,cspd,wdir,hd,wspd,tdir,tspd, + wmis) real crse(num), cspd(num), wdir(num), hd(num), wspd(num), + tdir(num), tspd(num), wmis(4) integer num print*,'\n-------------------------------------------------------- +----------------------------' print*,'\n MISSING DATA TABLE' print*,' ********************' print*,' index course sspeed windir shiphd appspd t +rudir truspd' do i=1,num if( (crse(i) .eq. wmis(1)) .or. (cspd(i) .eq. wmis(2)) .or. + (wdir(i) .eq. wmis(3)) .or. (wspd(i) .eq. wmis(4)) .or. + (hd(i) .eq. wmis(5)) ) then write(*,60) i, crse(i), cspd(i), wdir(i), hd(i), wspd(i), + tdir(i), tspd(i) 60 format(8X, i8, 7f8.1) endif enddo return end C ********************************************************************** c Subroutine: TRUERR.F c Purpose: List of where range tests fail and where values are c invalid. Also prints out number of records which are c complete, incomplete partially, incomplete entirely, and c where range tests fail. Accessed when selection #1, #2, c #3, or the default is chosen. subroutine truerr(num,crse,cspd,hd,wdir,wspd,wmis,nw,nwpm, + nwam,nwf) real crse(num), cspd(num), hd(num), wdir(num), wspd(num), + wmis(4) integer num, nw, nwpm, nwam, nwf character cv*22, wd*24, vs*22, ws*22, sh*15, str1*30, str2*10, + str3*25, str4*10 cv = 'Course value #' wd = 'Wind direction value #' vs = 'Vessel speed value #' ws = 'Wind speed value #' sh = 'Ship heading #' str1 = 'Truewinds range test failed. ' str2 = ' invalid.' str3 = 'Truewinds data test: ' str4 = ' missing.' print*,'\n-------------------------------------------------------- +----------------------------' print*,'\n TRUEWINDS ERRORS' print*,' ******************' do i=1, num if( ((crse(i) .lt. 0.0).or.(crse(i) .gt. 360.0)) .and. + (crse(i) .ne. wmis(1)) ) then write (*,70) str1, cv, i, str2 70 format(8X, A30, A14, i3, A9) endif if( (cspd(i) .lt. 0.0) .and. (cspd(i) .ne. wmis(2)) )then write(*,80) str1, vs, i, str2 80 format(8X, A30, A20, i3, A9) endif if( ((wdir(i) .lt. 0.0).or.(wdir(i) .gt. 360.0)) .and. + (wdir(i) .ne. wmis(3)) )then write (*,90) str1, wd, i, str2 90 format(8X, A30, A22, i3, A9) endif if( (wspd(i) .lt. 0.0) .and. (wspd(i) .ne. wmis(4)) ) then write(*,100) str1, ws, i, str2 100 format(8X, A30, A18, i3, A9) endif if( ((hd(i) .lt. 0.0).or.(hd(i) .gt. 360.0)) .and. + (hd(i) .ne. wmis(5)) )then write(*,110) str1, sh, i, str2 110 format(8X, A30, A14, i3, A9) endif enddo print *, ' ' do i=1, num if(crse(i) .eq. wmis(1)) then write(*,120) str3, cv, i, str4 120 format(8X, A22, A14, i3, A9) endif if(cspd(i) .eq. wmis(2)) then write(*,130) str3, vs, i, str4 130 format(8X, A22, A20, i3, A9) endif if(wdir(i) .eq. wmis(3)) then write(*,140) str3, wd, i, str4 140 format(8X, A22, A22, i3, A9) endif if(wspd(i) .eq. wmis(4)) then write(*,150) str3, ws, i, str4 150 format(8X, A22, A18, i3, A9) endif if(hd(i) .eq. wmis(5)) then write(*,160) str3, sh, i, str4 160 format(8X, A22, A14, I3, A9) endif enddo print*,'\n-------------------------------------------------------- +----------------------------' print*,'\n DATA REVIEW' print*,' *************' write(*,170) nw 170 format(26X, 'no data missing =',i4) write(*,180) nwpm 180 format(21X, 'part of data missing =',i4) write(*,190) nwam 190 format(25X, 'all data missing =',i4) write(*,200) nwf 200 format(23X, 'failed range tests =',i4) return end