program angles c c*********************************************************************c c* *c c* Written by Tim J. Cowan *c c* *c c*********************************************************************c c* c* Description: Converts Euler angles (roll, pitch, & yaw) to STARS c* CFD orientation angles (alpha & beta) for a specific c* free stream velocity vector in global coordinates. c* c* References: Peiro, Peraire, & Morgan, FELISA System Reference c* Manual, Part 2 - User Manual, December 1993: pg. 30-31. c* c* Nelson, Flight Stability and Automatic Control, 2nd c* Edition, McGraw-Hill, 1998: pg. 101-103. c* c* Comments: Watch-out for the miss-print in the Euler transformation c* matrix defined by Nelson! B(3,1) = sin(phi)*cos(theta) c* not sin(theta)*cos(theta). c* c* A velocity vector parallel to the x-axis is really the c* only choice that makes sense, i.e { 1, 0, 0 } or c* { -1, 0, 0 }. Otherwise the Euler angles aren't really c* appropriate. If you defined your geometry such that c* your free stream velocity vector is parallel to some c* other axis, you will need to think about your definition c* of the Euler angles for this to make sense. c* c*********************************************************************c c* c* Input : from user c* c* Output : to screen c* c*********************************************************************c c* c* Variable definitions c* c* name => description (default value) {additional information} c* c* alpha => STARS orientation angle #1 in degrees c* beta => STARS orientation angle #2 in degrees c* alp => STARS orientation angle #1 in radians c* bet => STARS orientation angle #2 in radians c* c* aeff => effective angle of attack in degrees c* beff => effective side slip angle in degrees c* c* roll => Euler angle #1 in degrees c* pitch => Euler angle #2 in degrees c* yaw => Euler angle #3 in degrees c* phi => Euler angle #1 in radians c* tet => Euler angle #2 in radians c* psi => Euler angle #3 in radians c* c* u,v,w => three components of free stream velocity vector c* u1,v1, c* & w1 => three components of rotated free stream velocity vector c* c* BT(9) => Euler transformation matrix c* BTI(9)=> inverse of Euler transformation matrix c* c*********************************************************************c c implicit real*8(A-H,O-Z) dimension BT(9), BTI(9) c c---------------------------------------------------------------------c c *** Problem definition : c---------------------------------------------------------------------c c pi = 4.0d0*DATAN(1.0d0) write(*,'(5(/),'' *** STARS Angle Calculator v2.1 ***'',3(/))') c c---------------------------------------------------------------------c c *** Get angles to convert from user : c---------------------------------------------------------------------c c write(*,'(''Enter a roll angle (phi) in degrees : '',$)') read(*,*) roll write(*,'(''Enter a pitch angle (theta) in degrees : '',$)') read(*,*) pitch write(*,'(''Enter a yaw angle (psi) in degrees : '',$)') read(*,*) yaw phi = roll*pi/180.0d0 tet = pitch*pi/180.0d0 psi = yaw*pi/180.0d0 c c---------------------------------------------------------------------c c *** Make sure the magnitude of the velocity vector is one : c---------------------------------------------------------------------c c u = 1.0d0 v = 0.0d0 w = 0.0d0 amag = DSQRT( u*u + v*v + w*w ) if ( amag .EQ. 0.0d0 ) then write(*,'(/,''Error: Your velocity vector is zero!'')') stop 'invalid velocity' endif u = u/amag v = v/amag w = w/amag c c---------------------------------------------------------------------c c *** Do the conversion : c---------------------------------------------------------------------c c write(*,'(2(/))') c *** assemble Euler angle transformation matrix & inverse c cp = dcos( phi ) sp = dsin( phi ) cq = dcos( tet ) sq = dsin( tet ) cr = dcos( psi ) sr = dsin( psi ) BT(1) = cq*cr BT(2) = sp*sq*cr - cp*sr BT(3) = cp*sq*cr + sp*sr BT(4) = cq*sr BT(5) = sp*sq*sr + cp*cr BT(6) = cp*sq*sr - sp*cr BT(7) = -sq BT(8) = sp*cq BT(9) = cp*cq BTI(1) = BT(5)*BT(9) - BT(6)*BT(8) BTI(2) = BT(3)*BT(8) - BT(2)*BT(9) BTI(3) = BT(2)*BT(6) - BT(3)*BT(5) BTI(4) = BT(6)*BT(7) - BT(4)*BT(9) BTI(5) = BT(1)*BT(9) - BT(7)*BT(3) BTI(6) = BT(3)*BT(4) - BT(1)*BT(6) BTI(7) = BT(4)*BT(8) - BT(5)*BT(7) BTI(8) = BT(2)*BT(7) - BT(1)*BT(8) BTI(9) = BT(1)*BT(5) - BT(2)*BT(4) c *** transform the free stream velocity vector c u1 = BTI(1)*u + BTI(2)*v + BTI(3)*w v1 = BTI(4)*u + BTI(5)*v + BTI(6)*w w1 = BTI(7)*u + BTI(8)*v + BTI(9)*w c *** calculate STARS angles c alp = DASIN( w1 ) bet = DASIN( v1/DCOS( alp ) ) alpha = alp*180.0d0/pi beta = bet*180.0d0/pi aeff = DATAN( w1/u1 )*180.0d0/pi beff = DATAN( v1/u1 )*180.0d0/pi c *** output results to screen c write(*,'(''STARS Orientation Angles : '')') write(*,'(''------------------------------------------'')') write(*,'(A,F12.6,A)') 'Alpha = ', alpha, ' deg.' write(*,'(A,F12.6,A)') 'Beta = ', beta, ' deg.' write(*,'(/)') write(*,'(''Euler Orientation Angles : '')') write(*,'(''------------------------------------------'')') write(*,'(A,F12.6,A)') 'Roll = ', roll, ' deg.' write(*,'(A,F12.6,A)') 'Pitch = ', pitch, ' deg.' write(*,'(A,F12.6,A)') 'Yaw = ', yaw, ' deg.' write(*,'(/)') write(*,'(''Aerodynamic Angles : '')') write(*,'(''------------------------------------------'')') write(*,'(A,F12.6,A)') 'Eff. angle of attack = ', aeff, ' deg.' write(*,'(A,F12.6,A)') 'Eff. side slip angle = ', beff, ' deg.' c c---------------------------------------------------------------------c c *** Clean up : c---------------------------------------------------------------------c c c c---------------------------------------------------------------------c c *** Successful completion of program : c---------------------------------------------------------------------c c 999 write (*,*) ' ' stop ' OK!!' c c---------------------------------------------------------------------c c *** Unsuccessful error messages : c---------------------------------------------------------------------c c end c*********************************************************************c c*********************************************************************c