1C*PGVECT -- vector map of a 2D data array, with blanking 2C%void cpgvect(const float *a, const float *b, int idim, int jdim, \ 3C% int i1, int i2, int j1, int j2, float c, int nc, \ 4C% const float *tr, float blank); 5C+ 6 SUBROUTINE PGVECT (A, B, IDIM, JDIM, I1, I2, J1, J2, C, NC, TR, 7 1 BLANK) 8 INTEGER IDIM, JDIM, I1, I2, J1, J2, NC 9 REAL A(IDIM,JDIM), B(IDIM, JDIM), TR(6), BLANK, C 10C 11C Draw a vector map of two arrays. This routine is similar to 12C PGCONB in that array elements that have the "magic value" defined by 13C the argument BLANK are ignored, making gaps in the vector map. The 14C routine may be useful for data measured on most but not all of the 15C points of a grid. Vectors are displayed as arrows; the style of the 16C arrowhead can be set with routine PGSAH, and the the size of the 17C arrowhead is determined by the current character size, set by PGSCH. 18C 19C Arguments: 20C A (input) : horizontal component data array. 21C B (input) : vertical component data array. 22C IDIM (input) : first dimension of A and B. 23C JDIM (input) : second dimension of A and B. 24C I1,I2 (input) : range of first index to be mapped (inclusive). 25C J1,J2 (input) : range of second index to be mapped (inclusive). 26C C (input) : scale factor for vector lengths, if 0.0, C will be 27C set so that the longest vector is equal to the 28C smaller of TR(2)+TR(3) and TR(5)+TR(6). 29C NC (input) : vector positioning code. 30C <0 vector head positioned on coordinates 31C >0 vector base positioned on coordinates 32C =0 vector centered on the coordinates 33C TR (input) : array defining a transformation between the I,J 34C grid of the array and the world coordinates. The 35C world coordinates of the array point A(I,J) are 36C given by: 37C X = TR(1) + TR(2)*I + TR(3)*J 38C Y = TR(4) + TR(5)*I + TR(6)*J 39C Usually TR(3) and TR(5) are zero - unless the 40C coordinate transformation involves a rotation 41C or shear. 42C BLANK (input) : elements of arrays A or B that are exactly equal to 43C this value are ignored (blanked). 44C-- 45C 4-Sep-1992: derived from PGCONB [J. Crane]. 46C 26-Nov-1992: revised to use PGARRO [TJP]. 47C 25-Mar-1994: correct error for NC not =0 [G. Gonczi]. 48C 5-Oct-1996: correct error in computing max vector length [TJP; 49C thanks to David Singleton]. 50C----------------------------------------------------------------------- 51 INTEGER I, J 52 REAL X, Y, X1, Y1, X2, Y2 53 REAL CC 54 INTRINSIC SQRT, MAX, MIN 55C 56C Define grid to world transformation 57C 58 X(I,J) = TR(1) + TR(2)*I + TR(3)*J 59 Y(I,J) = TR(4) + TR(5)*I + TR(6)*J 60C 61C Check arguments. 62C 63 IF (I1.LT.1 .OR. I2.GT.IDIM .OR. I1.GE.I2 .OR. 64 1 J1.LT.1 .OR. J2.GT.JDIM .OR. J1.GE.J2) THEN 65C CALL GRWARN('PGVECT: invalid range I1:I2, J1:J2') 66 RETURN 67 END IF 68C 69C Check for scale factor C. 70C 71 CC = C 72 IF (CC.EQ.0.0) THEN 73 DO 20 J=J1,J2 74 DO 10 I=I1,I2 75 IF (A(I,J).NE.BLANK .AND. B(I,J).NE.BLANK) 76 1 CC = MAX(CC,SQRT(A(I,J)**2+B(I,J)**2)) 77 10 CONTINUE 78 20 CONTINUE 79 IF (CC.EQ.0.0) RETURN 80 CC = SQRT(MIN(TR(2)**2+TR(3)**2,TR(5)**2+TR(6)**2))/CC 81 END IF 82C 83 CALL PGBBUF 84C 85 DO 40 J=J1,J2 86 DO 30 I=I1,I2 87C 88C Ignore vector if element of A and B are both equal to BLANK 89C 90 IF (.NOT.(A(I,J).EQ.BLANK .AND. B(I,J).EQ.BLANK)) THEN 91 92C 93C Define the vector starting and end points according to NC. 94C 95 IF (NC.LT.0) THEN 96 X2 = X(I,J) 97 Y2 = Y(I,J) 98 X1 = X2 - A(I,J)*CC 99 Y1 = Y2 - B(I,J)*CC 100 ELSE IF (NC.EQ.0) THEN 101 X2 = X(I,J) + 0.5*A(I,J)*CC 102 Y2 = Y(I,J) + 0.5*B(I,J)*CC 103 X1 = X2 - A(I,J)*CC 104 Y1 = Y2 - B(I,J)*CC 105 ELSE 106 X1 = X(I,J) 107 Y1 = Y(I,J) 108 X2 = X1 + A(I,J)*CC 109 Y2 = Y1 + B(I,J)*CC 110 END IF 111C 112C Draw vector. 113C 114 CALL PGARRO(X1, Y1, X2, Y2) 115 END IF 116 30 CONTINUE 117 40 CONTINUE 118C 119 CALL PGEBUF 120 END 121