WAVEWATCH III  beta 0.0.1
Data Types | Functions/Subroutines
w3cspcmd Module Reference

Convert spectra to new discrete spectral grid. More...

Data Types

type  case
 

Functions/Subroutines

subroutine w3cspc (SP1, NFR1, NTH1, XF1, FR1, TH1, SP2, NFR2, NTH2, XF2, FR2, TH2, NSP, NDST, NDSE, FTL)
 Convert a set of spectra to a new spectral grid. More...
 

Detailed Description

Convert spectra to new discrete spectral grid.

Author
H. L. Tolman
Date
01-Nov-2012
Copyright
Copyright 2009-2022 National Weather Service (NWS), National Oceanic and Atmospheric Administration. All rights reserved. WAVEWATCH III is a trademark of the NWS. No unauthorized use without permission.

Function/Subroutine Documentation

◆ w3cspc()

subroutine w3cspcmd::w3cspc ( real, dimension(nth1,nfr1,nsp), intent(in)  SP1,
integer, intent(in)  NFR1,
integer, intent(in)  NTH1,
real, intent(in)  XF1,
real, intent(in)  FR1,
real, intent(in)  TH1,
real, dimension(nth2,nfr2,nsp), intent(out)  SP2,
integer, intent(in)  NFR2,
integer, intent(in)  NTH2,
real, intent(in)  XF2,
real, intent(in)  FR2,
real, intent(in)  TH2,
integer, intent(in)  NSP,
integer, intent(in)  NDST,
integer, intent(in)  NDSE,
real, intent(in)  FTL 
)

Convert a set of spectra to a new spectral grid.

Conservative distribution of input energies over new grid.

Parameters
[in]SP1Input spectra.
[in]NFR1Input number of frequencies.
[in]NTH1Input number of directions.
[in]XF1Input frequency increment factor.
[in]FR1First input frequency.
[in]TH1First input direction.
[out]SP2Output spectra.
[in]NFR2Output number of frequencies.
[in]NTH2Output number of directions.
[in]XF2Output frequency increment factor.
[in]FR2First output frequency.
[in]TH2First output direction.
[in]NSPNumber of spectra.
[in]NDSTUnit number for test output.
[in]NDSEUnit number for error output.
[in]FTLFactor for tail description = XF2**N.
Author
H. L. Tolman
Date
01-Nov-2012

Definition at line 146 of file w3cspcmd.F90.

146  !/
147  !/ +-----------------------------------+
148  !/ | WAVEWATCH III NOAA/NCEP |
149  !/ | H. L. Tolman |
150  !/ | FORTRAN 90 |
151  !/ | Last update : 01-Nov-2012 !
152  !/ +-----------------------------------+
153  !/
154  !/ 19-Sep-2005 : Origination. ( version 3.08 )
155  !/ 01-Nov-2012 : code clean up (tab spaces, comments)( version 4.08 )
156  !/
157  ! 1. Purpose :
158  !
159  ! Convert a set of spectra to a new spectral grid.
160  !
161  ! 2. Method :
162  !
163  ! Conservative distribution of input energies over new grid.
164  !
165  ! 3. Parameters :
166  !
167  ! Parameter list
168  ! ----------------------------------------------------------------
169  ! SP1 R.A. I Input spectra.
170  ! NFR1 Int. I Input number of frequencies.
171  ! NTH1 Int. I Input number of directions.
172  ! XFR Real I Input frequency increment factor.
173  ! FR1 Real I First input frequency.
174  ! TH1 Real I First input direction.
175  ! SP2 R.A. O Output spectra.
176  ! NFR2, NTH2, XF2, FR2, TH2
177  ! ! Specral description for output spectra.
178  ! NSP Int. I Number of spectra.
179  ! NDST int. I Unit number for test output.
180  ! NDSE int. I Unit number for error output.
181  ! FTAIL Real I Factor for tail description = XF2**N
182  ! ----------------------------------------------------------------
183  !
184  ! 4. Subroutines used :
185  !
186  ! Name Type Module Description
187  ! ----------------------------------------------------------------
188  ! STRACE Sur. W3SERVMD Subroutine tracing.
189  ! EXTCDE Sur. Id program abort.
190  ! ----------------------------------------------------------------
191  !
192  ! 5. Called by :
193  !
194  ! Name Type Module Description
195  ! ----------------------------------------------------------------
196  ! W3IOBC Subr. W3IOBCMD Updating boundary conditions.
197  ! Subr Multi scale model bound. data input.
198  ! ----------------------------------------------------------------
199  !
200  ! 6. Error messages :
201  !
202  ! - Check on input parameters.
203  !
204  ! 7. Remarks :
205  !
206  ! - The inner loop of the actual redistribution is over the
207  ! individual spectra, optimizing this routine for large numbers
208  ! of conversions in a single call.
209  !
210  ! 8. Structure :
211  !
212  ! See source code.
213  !
214  ! 9. Switches :
215  !
216  ! !/S Enable subroutine tracing.
217  !
218  ! !/T Enable test output.
219  ! !/T1 Test output for searching in stored data.
220  ! !/T2 Test output for redistribution data.
221  !
222  ! 10. Source code :
223  !
224  !/ ------------------------------------------------------------------- /
225  USE constants
226  !
227  USE w3servmd, ONLY: extcde
228 #ifdef W3_S
229  USE w3servmd, ONLY: strace
230 #endif
231  !
232  IMPLICIT NONE
233  !/
234  !/ ------------------------------------------------------------------- /
235  !/ Parameter list
236  !/
237  INTEGER, INTENT(IN) :: NSP, NFR1, NTH1, NFR2, NTH2, NDST, NDSE
238  REAL, INTENT(IN) :: SP1(NTH1,NFR1,NSP), XF1, FR1, TH1, &
239  XF2, FR2, TH2, FTL
240  REAL, INTENT(OUT) :: SP2(NTH2,NFR2,NSP)
241  !/
242  !/ ------------------------------------------------------------------- /
243  !/ Local parameters
244  !/
245  INTEGER :: I, NRMAX, J, I1, L1, J1, I2, L2, J2, &
246  ISP
247 #ifdef W3_S
248  INTEGER, SAVE :: IENT = 0
249 #endif
250  REAL :: LOW, HGH, RLOW, RHGH, BLOW, BHGH, &
251  FRAC, AUX1, AUX2, R1, R2, FACT
252  LOGICAL :: FOUND
253  TYPE(CASE), POINTER :: CURRENT
254  !/
255  !/ ------------------------------------------------------------------- /
256  !/ Pointers for aliases
257  !/
258  INTEGER, POINTER :: IDTH(:,:), IDFR(:,:), NFR2T
259  REAL, POINTER :: DTH1, DTH2, RDTH(:,:), FRQ1(:), &
260  FRQ2(:), XDF1, XDF2, RDFR(:,:)
261  !/
262 #ifdef W3_S
263  CALL strace (ient, 'W3CSPC')
264 #endif
265  !
266  ! -------------------------------------------------------------------- /
267  ! 0. Initializations
268  ! 0.a Check input
269  !
270  IF ( nfr1.LT.3 .OR. nth1.LT.4 .OR. xf1.LE.1. .OR. fr1.LE.0. .OR.&
271  nfr2.LT.3 .OR. nth2.LT.4 .OR. xf2.LE.1. .OR. fr2.LE.0. ) THEN
272  WRITE (ndse,900) nfr1, nth1, xf1, fr1, nfr2, nth2, xf2, fr2
273  CALL extcde ( 1 )
274  END IF
275  !
276  IF ( nsp .LT. 0 ) THEN
277  WRITE (ndse,901)
278  CALL extcde ( 2 )
279  END IF
280  !
281  IF ( nsp .EQ. 0 ) THEN
282  WRITE (ndse,902)
283  RETURN
284  END IF
285  !
286  ! 0.b Test output
287  !
288 #ifdef W3_T
289  WRITE (ndst,9000) nsp, nfr1, nth1, xf1, fr1, th1*rade, &
290  nfr2, nth2, xf2, fr2, th2*rade, ftl
291 #endif
292  !
293  ! -------------------------------------------------------------------- /
294  ! 1. Search stored interpolation data for match
295  !
296  found = .false.
297  !
298  DO i=1, ncases
299  !
300  IF ( i .EQ. 1 ) THEN
301  current => idata
302  ELSE
303  current => current%NEXT
304  END IF
305  !
306 #ifdef W3_T1
307  WRITE (ndst,9010) i, current%NFR1, current%NTH1, &
308  current%XF1, current%FR1, current%TH1*rade, &
309  current%NFR2, current%NTH2, &
310  current%XF2, current%FR2, current%TH2*rade
311 #endif
312  !
313  found = current%NFR1.EQ.nfr1 .AND. current%NFR2.EQ.nfr2 .AND. &
314  current%NTH1.EQ.nth1 .AND. current%NTH2.EQ.nth2 .AND. &
315  current%XF1 .EQ.xf1 .AND. current%XF2 .EQ.xf2 .AND. &
316  current%FR1 .EQ.fr1 .AND. current%FR2 .EQ.fr2 .AND. &
317  current%TH1 .EQ.th1 .AND. current%TH2 .EQ.th2
318  IF ( found ) EXIT
319  !
320  END DO
321  !
322  ! -------------------------------------------------------------------- /
323  ! 2. Link or compute interpolation data
324  ! 2.a Link
325  !
326  IF ( found ) THEN
327  !
328 #ifdef W3_T
329  WRITE (ndst,9020) i
330 #endif
331  !
332  dth1 => current%DTH1
333  dth2 => current%DTH2
334  idth => current%IDTH
335  rdth => current%RDTH
336  !
337  frq1 => current%FRQ1
338  frq2 => current%FRQ2
339  xdf1 => current%XDF1
340  xdf2 => current%XDF2
341  nfr2t => current%NFR2T
342  idfr => current%IDFR
343  rdfr => current%RDFR
344  !
345  ! 2.b Compute
346  !
347  ELSE
348  !
349  ncases = ncases + 1
350 #ifdef W3_T
351  WRITE (ndst,9021) ncases
352 #endif
353  !
354  ! 2.b.1 Point and allocate as necessary
355  !
356  IF ( ncases .EQ. 1 ) THEN
357  ALLOCATE ( idata )
358  current => idata
359  ELSE
360  ALLOCATE ( current%NEXT )
361  current => current%NEXT
362  END IF
363  !
364  ! 2.b.2 Store test data
365  !
366  current%ICASE = ncases
367  current%NFR1 = nfr1
368  current%NTH1 = nth1
369  current%XF1 = xf1
370  current%FR1 = fr1
371  current%TH1 = th1
372  current%NFR2 = nfr2
373  current%NTH2 = nth2
374  current%XF2 = xf2
375  current%FR2 = fr2
376  current%TH2 = th2
377  !
378  ! 2.b.3 Directional redistribution data
379  !
380  dth1 => current%DTH1
381  dth1 = tpi / real(nth1)
382  dth2 => current%DTH2
383  dth2 = tpi / real(nth2)
384  !
385  IF ( dth1 .LE. dth2 ) THEN
386  nrmax = 2
387  ELSE
388  nrmax = 2 + int(dth1/dth2)
389  END IF
390  !
391  ALLOCATE (current%IDTH(0:nrmax,nth1),current%RDTH(nrmax,nth1))
392  idth => current%IDTH
393  rdth => current%RDTH
394  idth = 0
395  rdth = 0.
396  !
397  DO i=1, nth1
398  low = th1 + real(i-1)*dth1 - 0.5*dth1
399  hgh = low + dth1
400  rlow = 1. + (low-th2)/dth2
401  rhgh = 1. + (hgh-th2)/dth2
402  DO j=nint(rlow), nint(rlow)+nrmax-1
403  blow = th2 + real(j-1)*dth2 - 0.5*dth2
404  bhgh = blow + dth2
405  frac = (min(bhgh,hgh)-max(blow,low)) / (hgh-low)
406  IF ( frac .GT. 1.e-5 ) THEN
407  idth(0,i) = idth(0,i) + 1
408  idth(idth(0,i),i) = 1 + mod(j-1+nth2,nth2)
409  rdth(idth(0,i),i) = frac
410  END IF
411  END DO
412  END DO
413  !
414  ! 2.b.4 Frequency redistribution data
415  !
416  ALLOCATE ( current%FRQ1(nfr1), current%FRQ2(nfr2) )
417  frq1 => current%FRQ1
418  frq2 => current%FRQ2
419  !
420  frq1(1) = fr1
421  DO i=2, nfr1
422  frq1(i) = xf1 * frq1(i-1)
423  END DO
424  !
425  frq2(1) = fr2
426  DO i=2, nfr2
427  frq2(i) = xf2 * frq2(i-1)
428  END DO
429  !
430  xdf1 => current%XDF1
431  xdf1 = 0.5 * ( xf1 - 1./xf1 )
432  xdf2 => current%XDF2
433  xdf2 = 0.5 * ( xf2 - 1./xf2 )
434  !
435  IF ( xdf1 .LE. xdf2 ) THEN
436  nrmax = 2
437  ELSE
438  nrmax = 1
439  aux1 = xdf1
440  aux2 = xdf2
441  DO
442  nrmax = nrmax + 1
443  aux1 = aux1 - aux2
444  aux2 = aux2 / xf2
445  IF ( aux1 .LT. 0. ) EXIT
446  END DO
447  END IF
448  !
449  ALLOCATE (current%IDFR(0:nrmax,nfr1),current%RDFR(nrmax,nfr1))
450  idfr => current%IDFR
451  rdfr => current%RDFR
452  idfr = 0
453  rdfr = 0.
454  !
455  DO i=1, nfr1
456  IF ( i .EQ. 1 ) THEN
457  hgh = 0.5 * ( frq1(i) + frq1(i+1) )
458  low = hgh - xdf1*frq1(i)
459  ELSE
460  low = 0.5 * ( frq1(i) + frq1(i-1) )
461  hgh = low + xdf1*frq1(i)
462  END IF
463  DO j=1, nfr2
464  IF ( j .EQ. 1 ) THEN
465  bhgh = 0.5 * ( frq2(j) + frq2(j+1) )
466  blow = bhgh - xdf2*frq2(j)
467  ELSE
468  blow = 0.5 * ( frq2(j) + frq2(j-1) )
469  bhgh = blow + xdf2*frq2(j)
470  END IF
471  IF ( bhgh .LE. low ) cycle
472  IF ( blow .GE. hgh ) EXIT
473  frac = (min(bhgh,hgh)-max(blow,low)) / (hgh-low)
474  IF ( frac .LT. 1.e-5 ) cycle
475  idfr(0,i) = idfr(0,i) + 1
476  idfr(idfr(0,i),i) = j
477  rdfr(idfr(0,i),i) = frac
478  END DO
479  END DO
480  !
481  nfr2t => current%NFR2T
482  nfr2t = nfr2 + 1
483  DO j=nfr2, 1, -1
484  IF ( j .EQ. 1 ) THEN
485  bhgh = 0.5 * ( frq2(j) + frq2(j+1) )
486  ELSE
487  blow = 0.5 * ( frq2(j) + frq2(j-1) )
488  bhgh = blow + xdf2*frq2(j)
489  END IF
490  IF ( bhgh .GT. hgh ) THEN
491  nfr2t = j
492  ELSE
493  EXIT
494  END IF
495  END DO
496  !
497  END IF
498  !
499  ! 2.c Test output
500  !
501 #ifdef W3_T2
502  WRITE (ndst,9022)
503  DO i=1, nth1
504  WRITE (ndst,9024) i, idth(0,i), &
505  (idth(j,i),rdth(j,i),j=1,idth(0,i))
506  END DO
507  WRITE (ndst,9023) nfr2t
508  DO i=1, nfr1
509  WRITE (ndst,9024) i, idfr(0,i), &
510  (idfr(j,i),rdfr(j,i),j=1,idfr(0,i))
511  END DO
512 #endif
513  !
514  ! -------------------------------------------------------------------- /
515  ! 3. Convert
516  ! 3.a Discrete energies
517  !
518 #ifdef W3_T
519  WRITE (ndst,9030)
520 #endif
521  !
522  sp2 = 0.
523  !
524  DO i2=1, nfr1
525  DO l2=1, idfr(0,i2)
526  j2 = idfr(l2,i2)
527  r2 = rdfr(l2,i2)
528  DO i1=1,nth1
529  DO l1=1, idth( 0,i1)
530  j1 = idth(l1,i1)
531  r1 = rdth(l1,i1)
532  frac = r2 * frq1(i2) * xdf1 * r1 * dth1
533  sp2(j1,j2,:) = sp2(j1,j2,:) + frac * sp1(i1,i2,:)
534  END DO
535  END DO
536  END DO
537  END DO
538  !
539  ! 3.b Energy densities
540  !
541 #ifdef W3_T
542  WRITE (ndst,9031)
543 #endif
544  !
545  DO j2=1, nfr2
546  DO j1=1, nth2
547  fact = 1. / ( frq2(j2) * xdf2 * dth2 )
548  sp2(j1,j2,:) = fact * sp2(j1,j2,:)
549  END DO
550  END DO
551  !
552  ! 3.c Add the tail
553  !
554 #ifdef W3_T
555  WRITE (ndst,9032)
556 #endif
557  !
558  DO j2=nfr2t, nfr2
559  sp2(:,j2,:) = ftl * sp2(:,j2-1,:)
560  END DO
561  !
562  RETURN
563  !
564  ! Formats
565  !
566 900 FORMAT (/' *** ERROR W3CSPC: ILLEGAL INPUT PARAMETERS ***'/ &
567  ' INPUT : ',2i8,2f10.4/ &
568  ' OUTPUT : ',2i8,2f10.4)
569 901 FORMAT (/' *** ERROR W3CSPC: NEGATIVE NUMBER OF SPECTRA ***'/)
570 902 FORMAT (/' *** WARNING W3CSPC: NO SPECTRA ***'/)
571  !
572 #ifdef W3_T
573 9000 FORMAT ( ' TEST W3CSPC : NR. OF SPECTRA : ',i8/ &
574  ' INPUT SPECTRA : ',2i4,2f8.4,f6.1/ &
575  ' OUTPUT SPECTRA : ',2i4,2f8.4,f6.1/ &
576  ' TAIL FACTOR : ',f8.5)
577 #endif
578  !
579 #ifdef W3_T1
580 9010 FORMAT ( ' TEST W3CSPC : TEST INFO CASE : ',i8/ &
581  ' INPUT SPECTRA : ',2i4,2f8.4,f6.1/ &
582  ' OUTPUT SPECTRA : ',2i4,2f8.4,f6.1)
583 #endif
584  !
585 #ifdef W3_T
586 9020 FORMAT ( ' TEST W3CSPC : USING STORED DATA FOR CASE',i4)
587 9021 FORMAT ( ' TEST W3CSPC : COMPUTING DATA FOR CASE',i4)
588 #endif
589 #ifdef W3_T2
590 9022 FORMAT ( ' TEST W3CSPC : DIRECTIONAL DISTRIBUTION DATA')
591 9023 FORMAT ( ' TEST W3CSPC : FREQUENCY DISTRIBUTION DATA, ', &
592  'TAIL AT',i4)
593 9024 FORMAT ( ' ',i4,i4,' :',10(i4,f5.2) )
594 #endif
595  !
596 #ifdef W3_T
597 9030 FORMAT ( ' TEST W3CSPC : STARTING CONVERSION')
598 9031 FORMAT ( ' TEST W3CSPC : ENERGIES TO DENSITIES')
599 9032 FORMAT ( ' TEST W3CSPC : ADD TAIL')
600 #endif
601  !/
602  !/ End of W3CSPC ----------------------------------------------------- /
603  !/

References w3servmd::extcde(), constants::rade, w3servmd::strace(), and constants::tpi.

Referenced by w3bounc(), w3iobcmd::w3iobc(), wminiomd::wmiobg(), wminiomd::wmiobs(), wminiomd::wmioeg(), wminiomd::wmiohg(), and wmiopomd::wmiopo().

constants::rade
real, parameter rade
RADE Conversion factor from radians to degrees.
Definition: constants.F90:76
w3servmd
Definition: w3servmd.F90:3
constants::tpi
real, parameter tpi
TPI 2*Pi.
Definition: constants.F90:72
w3servmd::strace
subroutine strace(IENT, SNAME)
Definition: w3servmd.F90:148
constants
Define some much-used constants for global use (all defined as PARAMETER).
Definition: constants.F90:20
w3servmd::extcde
subroutine extcde(IEXIT, UNIT, MSG, FILE, LINE, COMM)
Definition: w3servmd.F90:736