-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathwwm_bdcons.F90
3068 lines (2935 loc) · 114 KB
/
wwm_bdcons.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
#include "wwm_functions.h"
#define DEBUG
#undef DEBUG
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE SPECTRAL_SHAPE(SPPAR,WALOC,LDEBUG,CALLFROM, OPTI)
USE DATAPOOL
IMPLICIT NONE
REAL(rkind), INTENT(OUT) :: WALOC(NUMSIG,NUMDIR)
REAL(rkind), INTENT(INOUT) :: SPPAR(8)
CHARACTER(LEN=*), INTENT(IN) :: CALLFROM
LOGICAL, INTENT(IN) :: LDEBUG, OPTI
IF (SPPAR(1) .lt. VERYSMALL) THEN
CALL KERNEL_SPECTRAL_SHAPE(SPPAR,WALOC,LDEBUG,CALLFROM)
ELSE
IF (OPTI) THEN
! Print *, 'Before call to OPTI_SPECTRAL_SHAPE'
CALL OPTI_SPECTRAL_SHAPE(SPPAR,WALOC,LDEBUG,CALLFROM)
! Print *, ' After call to OPTI_SPECTRAL_SHAPE'
ELSE
CALL KERNEL_SPECTRAL_SHAPE(SPPAR,WALOC,LDEBUG,CALLFROM)
END IF
END IF
END SUBROUTINE
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE COMPUTE_ESTIMATE_PER_DIR_SHAPE(SPPAR, WALOC, HS, TM, DM)
USE DATAPOOL
IMPLICIT NONE
REAL(rkind), INTENT(IN) :: SPPAR(8)
REAL(rkind), INTENT(IN) :: WALOC(NUMSIG,NUMDIR)
REAL(rkind), INTENT(OUT) :: HS, TM, DM
REAL(rkind) :: DEPLOC, CURTXYLOC(2)
REAL(rkind) :: WKLOC(NUMSIG)
REAL(rkind) :: FPP, TPP, CPP, WNPP, CGPP, KPP, LPP
REAL(rkind) :: PEAKDSPR, PEAKDM, DPEAK, TPPD, KPPD, CGPD, CPPD
REAL(rkind) :: TM01, TM02, TM10, KLM, WLM
REAL(rkind) :: ETOTS, ETOTC, DSPR
REAL(rkind) :: SPSIGLOC, WVN, WVC, WVK, WVCG
integer ISMAX, IS
ISMAX=NUMSIG
CURTXYLOC=ZERO
DO IS=1,NUMSIG
SPSIGLOC = SPSIG(IS)
CALL WAVEKCG(DEPLOC,SPSIGLOC,WVN,WVC,WVK,WVCG)
WKLOC(IS)=WVK
END DO
CALL MEAN_PARAMETER_LOC(WALOC,CURTXYLOC,DEPLOC,WKLOC,ISMAX,HS,TM01,TM02,TM10,KLM,WLM)
IF (SPPAR(5) .gt. 0) THEN
! Print *, 'Using PEAK parameters'
CALL PEAK_PARAMETER_LOC(WALOC,DEPLOC,ISMAX,FPP,TPP,CPP,WNPP,CGPP,KPP,LPP,PEAKDSPR,PEAKDM,DPEAK,TPPD,KPPD,CGPD,CPPD)
! Print *, ' PEAKDM=', PEAKDM
TM=TPP
DM=PEAKDM
ELSE
! Print *, 'Using MEAN parameters'
CALL MEAN_DIRECTION_AND_SPREAD_LOC(WALOC,ISMAX,ETOTS,ETOTC,DM,DSPR)
TM=TM01
END IF
END SUBROUTINE
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE OPTI_SPECTRAL_SHAPE(SPPAR,WALOC,LDEBUG,CALLFROM)
USE DATAPOOL
IMPLICIT NONE
REAL(rkind), INTENT(OUT) :: WALOC(NUMSIG,NUMDIR)
REAL(rkind), INTENT(INOUT) :: SPPAR(8)
CHARACTER(LEN=*), INTENT(IN) :: CALLFROM
LOGICAL, INTENT(IN) :: LDEBUG
REAL(rkind) :: HS, TM, DM, TheErr, DeltaPer, Tper
REAL(rkind) :: DiffAng, DEG, ADIR
REAL(rkind) :: SPPARwork1(8), SPPARwork2(8), SPPARwork(8)
integer :: iIter, nbIter, eSign, IS, ID
REAL(rkind) :: eSum
IF (ABS(SPPAR(5)) .eq. 3) THEN
CALL KERNEL_SPECTRAL_SHAPE(SPPAR,WALOC,LDEBUG,CALLFROM)
RETURN
END IF
SPPARwork1=SPPAR
SPPARwork2=SPPAR
Tper=SPPAR(2)
CALL KERNEL_SPECTRAL_SHAPE(SPPAR,WALOC,LDEBUG,CALLFROM)
CALL COMPUTE_ESTIMATE_PER_DIR_SHAPE(SPPAR, WALOC, HS, TM, DM)
DeltaPer=Tper - TM
IF (TM < Tper) THEN
eSign=1
ELSE
eSign=-1
END IF
SPPARwork=SPPAR
SPPARwork(2)=SPPAR(2) + DeltaPer
iIter=0
DO
iIter=iIter + 1
CALL KERNEL_SPECTRAL_SHAPE(SPPARwork,WALOC,LDEBUG,CALLFROM)
CALL COMPUTE_ESTIMATE_PER_DIR_SHAPE(SPPAR, WALOC, HS, TM, DM)
TheErr=(TM - Tper)*eSign
! Print *, 'iIter=', iIter, ' TheErr=', TheErr, ' DeltaPer=', DeltaPer
! Print *, ' eSign=', eSign, ' TM=', TM, ' Tper=', Tper
! Print *, 'SPPARwork(2)=', SPPARwork(2)
IF (TheErr > 0) THEN
EXIT
END IF
SPPARwork(2)=SPPARwork(2) + DeltaPer
END DO
IF (eSign .eq. 1) THEN
SPPARwork1=SPPAR
SPPARwork2=SPPARwork
ELSE
SPPARwork1=SPPARwork
SPPARwork2=SPPAR
END IF
nbIter=20
iIter=0
DO
! Print *, 'iIter=', iIter
SPPARwork=0.5_rkind*SPPARwork1 + 0.5_rkind*SPPARwork2
CALL KERNEL_SPECTRAL_SHAPE(SPPARwork,WALOC,LDEBUG,CALLFROM)
CALL COMPUTE_ESTIMATE_PER_DIR_SHAPE(SPPAR, WALOC, HS, TM, DM)
IF (TM > Tper) THEN
SPPARwork2=SPPARwork
ELSE
SPPARwork1=SPPARwork
END IF
iIter=iIter + 1
IF (iIter > nbIter) THEN
EXIT
END IF
END DO
CALL KERNEL_SPECTRAL_SHAPE(SPPARwork,WALOC,LDEBUG,CALLFROM)
CALL COMPUTE_ESTIMATE_PER_DIR_SHAPE(SPPAR, WALOC, HS, TM, DM)
SPPARwork(1)=SPPAR(1)*(SPPAR(1)/HS)
CALL KERNEL_SPECTRAL_SHAPE(SPPARwork,WALOC,LDEBUG,CALLFROM)
DO IS=1,NUMSIG
eSum=sum(WALOC(IS,:))
END DO
CALL DEG2NAUT (SPPAR(3), DEG, LNAUTIN)
ADIR = DEG * DEGRAD
DO ID=1,NUMDIR
eSum=sum(WALOC(:,ID))
DiffAng=(360.0_rkind/PI2)*(SPDIR(ID) - ADIR)
END DO
END SUBROUTINE
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE KERNEL_SPECTRAL_SHAPE(SPPAR,WALOC,LDEBUG,CALLFROM)
USE DATAPOOL
IMPLICIT NONE
REAL(rkind), INTENT(OUT) :: WALOC(NUMSIG,NUMDIR)
REAL(rkind), INTENT(INOUT) :: SPPAR(8)
CHARACTER(LEN=*), INTENT(IN) :: CALLFROM
LOGICAL, INTENT(IN) :: LDEBUG
INTEGER ID, IS, LSHAPE, ITPER, ISPEAK, ISIGMP
REAL(rkind) :: APSHAP, AUX1, AUX2, AUX3, AM0, AM1, AS2, AS3, ETOT, VEC2DEG
REAL(rkind) :: COEFF, SYF , MPER, PKPER, DIFPER, EHFR
REAL(rkind) :: MS, DEG, ETOTS, ETOTC, FF, CPSHAP, PPSHAP, DM, EAD, DS
REAL(rkind) :: RA, SALPHA, SF, SF4, SF5, FPK, FPK4, EFTAIL, CDIR
REAL(rkind) :: GAMMA_FUNC, DSPR, AACOS, ADIR, ETAIL_ARR, ATAIL_ARR, PTAIL_ARR
REAL(rkind) :: OMEG, EFTOT, ETOTT, CTOT, TM1, TPEAK
LOGICAL LOGPM, LINCOUT
! SPPARM(1), WBHS: Hs, sign. wave height
! SPPARM(2), WBTP: Wave period given by the user (either peak or mean)
! SPPARM(3), WBDM: average direction
! SPPARM(4), WBDS: directional spread
! SPPARM(5), WBSS: spectral shape (1-4), (1 - Pierson-Moskowitz, 2 - JONSWAP, 3 - BIN, 4 - Gauss) peak (+) or mean frequency (-)
! SPPARM(6), WBDSMS: directional spreading in degree (2) or exponent (1)
! SPPARM(7), WBGAUSS: gaussian width for the gauss spectrum 0.1
! SPPARM(8), WBPKEN: peak enhancement factor for the JONSWAP spectra 3.3
IF (LDEBUG) THEN
WRITE(DBG%FHNDL,*) 'HS PER DIR DPSR SHAPE DEGEXP GAUSS PEAK'
WRITE(DBG%FHNDL,'(8F10.4)') SPPAR(8)
ENDIF
ETOT = 0.
EFTOT = 0.
WALOC = 0.
IF (SPPAR(1) .LT. THR .OR. SPPAR(2) .LT. THR .OR. SPPAR(4) .LT. THR) THEN
WALOC = 0.
RETURN
END IF
IF (SPPAR(5) .LT. 0) THEN
LSHAPE = -INT(SPPAR(5))
LOGPM = .FALSE.
ELSE
LSHAPE = INT(SPPAR(5))
LOGPM = .TRUE.
ENDIF
!
PKPER = SPPAR(2)
ITPER = 0
IF (LSHAPE.EQ.3) THEN
! select bin closest to given period
DIFPER = 1.E10
DO IS = 1, NUMSIG
IF (ABS(PKPER - PI2/SPSIG(IS)) .LT. DIFPER) THEN
ISPEAK = IS
DIFPER = ABS(PKPER - PI2/SPSIG(IS))
END IF
ENDDO
ENDIF
!
100 FPK = (1./PKPER)
FPK4 = FPK**4
IF (LSHAPE.EQ.1) THEN
SALPHA = SPPAR(1)**2*FPK4 * 5./16.
ELSE IF (LSHAPE.EQ.2) THEN
SALPHA = (SPPAR(1)**2 * FPK4) / ((0.06533*(SPPAR(8)**0.8015)+0.13467)*16.)
ELSE IF (LSHAPE.EQ.4) THEN
AUX1 = SPPAR(1)**2 / ( 16.* SQRT (PI2) * SPPAR(7))
AUX3 = 2._rkind * SPPAR(7)**2
ENDIF
!
DO IS = 1, NUMSIG
!
IF (LSHAPE.EQ.1) THEN
SF = SPSIG(IS) / PI2
SF4 = SF**4
SF5 = SF**5
RA = (SALPHA/SF5)*EXP(-(5.*FPK4)/(4.*SF4))/(PI2*SPSIG(IS))
WALOC(IS,NUMDIR) = RA
ELSE IF (LSHAPE.EQ.2) THEN
SF = SPSIG(IS)/(PI2)
SF4 = SF**4
SF5 = SF**5
CPSHAP = 1.25_rkind * FPK4 / SF4
IF (CPSHAP.GT.10._rkind) THEN
RA = 0._rkind
ELSE
RA = (SALPHA/SF5) * EXP(-CPSHAP)
ENDIF
IF (SF .LT. FPK) THEN
COEFF = 0.07_rkind
ELSE
COEFF = 0.09_rkind
ENDIF
IF (COEFF*FPK .GT. SMALL) THEN
APSHAP = 0.5_rkind * ((SF-FPK) / (COEFF*FPK))**2
ELSE
APSHAP = ZERO
ENDIF
IF (APSHAP .GT. 10._rkind) THEN
SYF = 1.
ELSE
PPSHAP = EXP(-APSHAP)
SYF = SPPAR(8)**PPSHAP
ENDIF
RA = SYF*RA/(SPSIG(IS)*PI2)
WALOC(IS,NUMDIR) = RA
IF (LDEBUG) WRITE(DBG%FHNDL,*) 'IS LOOP', IS, SF, FPK, SYF, RA
!
ELSE IF (LSHAPE .EQ. 3) THEN
IF (IS.EQ.ISPEAK) THEN
ISBIN = ISPEAK
WALOC(IS,NUMDIR) = ( SPPAR(1)**2 ) / ( 16. * SIGPOW(IS,2) * FRINTF )
ELSE
WALOC(IS,NUMDIR) = 0.
END IF
ELSE IF (LSHAPE .EQ. 4) THEN
AUX2 = ( SPSIG(IS) - ( PI2 / PKPER ) )**2
RA = AUX1 * EXP ( -1. * AUX2 / AUX3 ) / SPSIG(IS)
WALOC(IS,NUMDIR) = RA
ELSE
CALL WWM_ABORT('Wrong type for frequency shape 1 - 4')
ENDIF
END DO
MPER = 0.
IF (.NOT.LOGPM.AND.ITPER.LT.100) THEN
ITPER = ITPER + 1
! calculate average frequency
AM0 = 0.
AM1 = 0.
DO IS = 1, NUMSIG
AS2 = WALOC(IS,NUMDIR) * SIGPOW(IS,2)
AS3 = AS2 * SPSIG(IS)
AM0 = AM0 + AS2
AM1 = AM1 + AS3
ENDDO
! contribution of tail to total energy density
PTAIL_ARR = TAIL_ARR(1) - 1.
ATAIL_ARR = 1. / (PTAIL_ARR * (1. + PTAIL_ARR * (FRINTH-1.)))
AM0 = AM0 * FRINTF + ATAIL_ARR * AS2
PTAIL_ARR = TAIL_ARR(1) - 2.
ETAIL_ARR = 1. / (PTAIL_ARR * (1. + PTAIL_ARR * (FRINTH-1.)))
AM1 = AM1 * FRINTF + ETAIL_ARR * AS3
! Mean period:
IF ( AM1.GT.THR) THEN
MPER = PI2 * AM0 / AM1
ELSE
MPER = ZERO
ENDIF
! write(*,'(I10,4F15.8)') ITPER, MPER, &
! & ABS(MPER-SPPAR(2)) .GT. 0.01*SPPAR(2), &
! & (SPPAR(2) / MPER) * PKPER
IF (ABS(MPER-SPPAR(2)) .GT. 0.01*SPPAR(2) .AND. MPER .GT. THR) THEN
! modification suggested by Mauro Sclavo
PKPER = (SPPAR(2)/MPER) * PKPER
GOTO 100
ELSE
PKPER = ZERO
ENDIF
ELSE IF (ITPER.GE.100) THEN
WRITE(STAT%FHNDL,*) 'No convergence calculating the spectrum'
FLUSH(STAT%FHNDL)
ENDIF
CALL DEG2NAUT (SPPAR(3), DEG, LNAUTIN)
ADIR = DEG * DEGRAD
IF (INT(SPPAR(6)) .EQ. 1) THEN
DSPR = PI * SPPAR(4) / 180._rkind
MS = MAX (DSPR**(-2) - TWO, 1._rkind)
ELSE
MS = SPPAR(4)
ENDIF
IF (MS.LT.12._rkind) THEN
CTOT = (2._rkind**MS)*(GAMMA_FUNC(0.5_rkind*MS+1._rkind))**2/(PI*GAMMA_FUNC(MS+1._rkind))
ELSE
CTOT = SQRT(0.5_rkind*MS/PI)/(1._rkind-0.25_rkind/MS)
ENDIF
LINCOUT = .FALSE.
DO ID = 1, NUMDIR
AACOS = COS(SPDIR(ID) - ADIR)
!write(*,*) aacos, spdir(id), adir
IF (AACOS .GT. 0._rkind) THEN
CDIR = CTOT * MAX (AACOS**MS, THR)
IF (.NOT.LCIRD) THEN
IF (AACOS .GE. COS(DDIR)) THEN
LINCOUT = .TRUE.
!WRITE(*,*) 'ERROR', AACOS, COS(DDIR)
!STOP 'AVERAGE DIRECTION IS OUT OF SECTOR'
ENDIF
ENDIF
ELSE
CDIR = 0._rkind
ENDIF
DO IS = 1, NUMSIG
WALOC(IS,ID) = CDIR * WALOC(IS,NUMDIR)
!write(*,'(2I10,2F15.8)') is, id, cdir, WALOC(IS,NUMDIR)
ENDDO
ENDDO
!
! Finished creating
!
! Integral Parameters of the Input Spectra ...
! AR: 2DO: Optimize the output routine for this task ...
!
IF (LDEBUG) THEN
ETOT = 0.
ETOTC = 0.
ETOTS = 0.
EFTAIL = 1.0 / (TAIL_ARR(1)-1.0)
IF (NUMSIG .GE. 2) THEN
DO ID = 1, NUMDIR
DO IS = 2, NUMSIG
DS = SPSIG(IS) - SPSIG(IS-1)
EAD = 0.5*(SPSIG(IS)*WALOC(IS,ID)+SPSIG(IS-1)*WALOC(IS-1,ID))*DS*DDIR
ETOT = ETOT + EAD
ETOTC = ETOTC + EAD * COS(SPDIR(ID))
ETOTS = ETOTS + EAD * SIN(SPDIR(ID))
END DO
IF (NUMSIG > 3) THEN
EHFR = WALOC(NUMSIG,ID) * SPSIG(NUMSIG)
ETOT = ETOT + DDIR * EHFR * SPSIG(NUMSIG) * EFTAIL
ENDIF
END DO
ELSE
DS = SGHIGH - SGLOW
DO ID = 1, NUMDIR
EAD = WALOC(1,ID) * DS * DDIR
ETOT = ETOT + EAD
ETOTC = ETOTC + EAD * COS(SPDIR(ID))
ETOTS = ETOTS + EAD * SIN(SPDIR(ID))
END DO
END IF
IF (ETOT > THR ) THEN
DM = VEC2DEG (ETOTC, ETOTS)
CALL DEG2NAUT(DM,DEG,LNAUTOUT)
DM = DEG
FF = MIN (ONE, SQRT(ETOTC*ETOTC+ETOTS*ETOTS)/ETOT)
DSPR = SQRT(2.-2.*FF) * 180./PI
ELSE
DM = 0.
DSPR = 0.
END IF
ETOTT = 0.0
EFTOT = 0.0
EFTAIL = TAIL_ARR(3)
DO ID = 1, NUMDIR
DO IS = 1, NUMSIG
OMEG = SPSIG(IS)
EAD = FRINTF * SIGPOW(IS,2) * WALOC(IS,ID)
ETOTT = ETOTT + EAD
EFTOT = EFTOT + EAD * OMEG
END DO
END DO
IF (EFTOT > VERYSMALL) THEN
TM1 = PI2*ETOTT/EFTOT
ELSE
TM1 = 0.0
END IF
ETOTT = 0.0
ISIGMP = -1
DO IS = 1, NUMSIG
EAD = 0.0
DO ID = 1, NUMDIR
EAD = EAD + SPSIG(IS)*WALOC(IS,ID)*DDIR
ENDDO
IF (EAD > ETOTT) THEN
ETOTT = EAD
ISIGMP = IS
END IF
END DO
IF (ISIGMP > 0) THEN
TPEAK = 1./(SPSIG(ISIGMP)/PI2)
ELSE
TPEAK = 0.0
END IF
IF (PrintLOG) THEN
WRITE(STAT%FHNDL,'("+TRACE...",A)') 'GIVEN BOUNDARY SPECTRA AND RECALCULATED WAVE PARAMETERS'
WRITE(STAT%FHNDL,'("+TRACE...",A)') 'THE DIFFERENCE IS MOSTLY DUE TO COARSE RESOLUTION IN SPECTRAL SPACE'
WRITE(STAT%FHNDL,*) 'GIVEN ', 'HS =', SPPAR(1)
WRITE(STAT%FHNDL,*) 'GIVEN ', 'DM =', SPPAR(3)
WRITE(STAT%FHNDL,*) 'GIVEN ', 'DSPR =', SPPAR(4)
WRITE(STAT%FHNDL,*) 'GIVEN ', 'TM or TP', SPPAR(2)
WRITE(STAT%FHNDL,*) 'SIMUL ', 'HS =', 4. * SQRT(ETOT)
WRITE(STAT%FHNDL,*) 'SIMUL ', 'DM =', DM
WRITE(STAT%FHNDL,*) 'SIMUL ', 'DSPR =', DSPR
WRITE(STAT%FHNDL,*) 'SIMUL ', 'TM=', TM1, 'TPEAK=', TPEAK
WRITE(STAT%FHNDL,*) 'TOT AC =', SUM(WALOC)
WRITE(STAT%FHNDL,*) SPPAR
FLUSH(STAT%FHNDL)
END IF
END IF
END SUBROUTINE
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE SPECTRUM_INT(WALOC)
USE DATAPOOL
IMPLICIT NONE
REAL(rkind), INTENT(OUT) :: WALOC(NUMSIG,NUMDIR)
REAL(rkind) :: MS(NUMSIG), MS1, ADIR1, DS, EAD
REAL(rkind) :: INSPF(WBNUMSIG)
REAL(rkind) :: INSPE(WBNUMSIG)
REAL(rkind) :: INDIR(WBNUMSIG)
REAL(rkind) :: INSPRD(WBNUMSIG)
REAL(rkind) :: INMS(WBNUMSIG)
REAL(rkind) :: SPCDIR(NUMSIG)
INTEGER :: IS, IS2, ID
REAL(rkind) :: CTOT(NUMSIG), CDIRT, CDIR(NUMDIR), CTOT1, CDIR1
REAL(rkind) :: DDACOS, DEG, DX, DIFFDX, YINTER
REAL(rkind) :: GAMMA_FUNC, ETOT, TM2
REAL(rkind) :: EFTOT, TM1, OMEG, PTAIL_ARR, OMEG2
REAL(rkind) :: RA, ETAIL, EFTAIL
REAL(rkind) :: THD(NUMDIR)
REAL(rkind) :: DTHD, RTH0
MS1 = WBDS
ADIR1 = WBDM
PSHAP(1) = 3.3
PSHAP(2) = 0.07
PSHAP(3) = 0.09
PSHAP(4) = 0.01
PSHAP(5) = 2.0
PSHAP(6) = 0.01
!
! 2do Convert ... set the correct units ...
!
DO IS = 1, WBNUMSIG
INSPF(IS) = SFRQ(IS,1) * PI2
INDIR(IS) = SDIR(IS,1)
INSPRD(IS) = SPRD(IS,1) * DEGRAD
INSPE(IS) = SPEG(IS,1,1) / PI2
END DO
ETOT = 0.0
DO IS = 2, WBNUMSIG
DS = (INSPF(IS) - INSPF(IS-1))
EAD = 0.5*(INSPE(IS) + INSPE(IS-1))*DS
ETOT = ETOT + EAD
END DO
IF (PrintLOG) THEN
WRITE (STAT%FHNDL,*) 'HS - INPUTSPECTRA - 1', 4.0*SQRT(ETOT)
END IF
WALOC = 0.
CALL INTERLIN (WBNUMSIG, NUMSIG, INSPF, SPSIG, INSPE, WALOC(:,1))
DO IS = 1, NUMSIG
IF (SPSIG(IS) .GT. INSPF(WBNUMSIG)) THEN
IF (PrintLOG) THEN
WRITE(STAT%FHNDL,*) 'Discrete Frequency is bigger then measured set FRMAX =', INSPF(WBNUMSIG)/PI2
WRITE(STAT%FHNDL,*) 'Setting all Action above the max. measured freq. zero'
END IF
WALOC(IS,1) = 0.0
END IF
END DO
ETOT = 0.0
DO IS = 2, NUMSIG
DS = SPSIG(IS) - SPSIG(IS-1)
EAD = 0.5*(WALOC(IS,1) + WALOC(IS-1,1))*DS
ETOT = ETOT + EAD
END DO
IF (PrintLOG) THEN
WRITE (STAT%FHNDL,*) 'HS - INPUTSPECTRA - INTERPOLATED', 4.0*SQRT(ETOT)
END IF
!
! Convert to Wave Action if nessasary
!
IF (LWBAC2EN) THEN
DO IS = 1, NUMSIG
WALOC(IS,1) = WALOC(IS,1) / SPSIG(IS)
END DO
END IF
ETOT = 0.0
DO IS = 2, NUMSIG
DS = SPSIG(IS) - SPSIG(IS-1)
EAD = 0.5 * (SPSIG(IS)*WALOC(IS,1)+SPSIG(IS-1)*WALOC(IS-1,1))*DS
ETOT = ETOT + EAD
END DO
IF (PrintLOG) THEN
WRITE(STAT%FHNDL,*) 'HS - INPUTSPECTRA - WAVE ACTION', 4.0*SQRT(ETOT)
END IF
!
! Convert from nautical to cartesian direction if nessasary and from deg2rad
!
DO IS = 1, WBNUMSIG
CALL DEG2NAUT (INDIR(IS), DEG, LNAUTIN)
INDIR(IS) = DEG
END DO
!
! Interpolate Directions in frequency space
!
DO IS = 1, NUMSIG
DO IS2 = 1, WBNUMSIG - 1
IF (SPSIG(IS) .GT. INSPF(IS2) .AND. SPSIG(IS) .LT. INSPF(IS2+1)) THEN
DX = INSPF(IS2+1) - INSPF(IS2)
DIFFDX = SPSIG(IS) - INSPF(IS2)
CALL INTERDIR( INDIR(IS2), INDIR(IS2+1), DX, DIFFDX, YINTER)
SPCDIR(IS) = YINTER
IF (SPSIG(IS) .GT. INSPF(WBNUMSIG) ) SPCDIR(IS) = 0.0
END IF
END DO
IF (SPSIG(IS) .GT. INSPF(WBNUMSIG) ) SPCDIR(IS) = 0.0
END DO
CALL INTERLIN (WBNUMSIG, NUMSIG, INSPF, SPSIG, INDIR, SPCDIR)
DO IS = 1, NUMSIG
IF ( SPSIG(IS) .GT. INSPF(WBNUMSIG) ) SPCDIR(IS) = 0.
END DO
DO IS = 1, NUMSIG
DEG = SPCDIR(IS) * DEGRAD
SPCDIR(IS) = DEG
END DO
IF (LINDSPRDEG) THEN
DO IS = 1, WBNUMSIG
INMS(IS) = MAX (INSPRD(IS)**(-2) - TWO, ONE)
END DO
ELSE
DO IS = 1, WBNUMSIG
INMS(IS) = INSPRD(IS)
END DO
END IF
!
! Interpolate MS in Frequency Space, if LCUBIC than Cubic Spline Interpolation is used
!
MS = 0.
CALL INTERLIN (WBNUMSIG, NUMSIG, INSPF, SPSIG, INMS, MS)
DO IS = 1, NUMSIG
IF ( SPSIG(IS) .GT. INSPF(WBNUMSIG) ) MS(IS) = MS(IS-1)
END DO
!
! Construction of a JONSWAP Spectra
!
IF (LPARMDIR) THEN
IF (MS1.GT.10.) THEN
CTOT1 = SQRT(MS1/(2.*PI)) * (1. + 0.25/MS1)
ELSE
CTOT1 = 2.**MS1 * (GAMMA_FUNC(1.+0.5*MS1))**2 / (PI * GAMMA_FUNC(1.+MS1))
ENDIF
DO IS = 1, NUMSIG
RA = WALOC(IS,1)
CDIRT = 0.
DO ID = 1, NUMDIR
DDACOS = COS(SPDIR(ID) - ADIR1)
IF (DDACOS .GT. 0.) THEN
CDIR1 = CTOT1 * MAX (DDACOS**MS1, THR)
ELSE
CDIR1 = 0.
ENDIF
WALOC(IS,ID) = RA * CDIR1
CDIRT = CDIRT + CDIR1 * DDIR
END DO
END DO
ELSE
DO IS = 1, NUMSIG
IF (MS(IS).GT.10.) THEN
CTOT(IS) = SQRT(MS(IS)/(2.*PI)) * (1. + 0.25/MS(IS))
ELSE
CTOT(IS) = 2.**MS(IS) * (GAMMA_FUNC(1.+0.5*MS(IS)))**2.0 / (PI * GAMMA_FUNC(1.+MS(IS)))
ENDIF
END DO
DO IS = 1, NUMSIG
RA = WALOC(IS,1)
CDIRT = 0.
DO ID = 1, NUMDIR
DDACOS = COS(SPDIR(ID) - SPCDIR(IS))
IF (DDACOS .GT. 0.) THEN
CDIR(ID) = CTOT(IS) * MAX (DDACOS**MS(IS), ZERO)
ELSE
CDIR(ID) = 0.
ENDIF
WALOC(IS,ID) = RA * CDIR(ID)
CDIRT = CDIRT + CDIR(ID) * DDIR
ENDDO
IF (100. - 1./CDIRT*100. .GT. 1.) THEN
WRITE (STAT%FHNDL,*) 100 - 1./CDIRT*100., 'ERROR BIGGER THAN 1% IN THE DIRECTIONAL DISTTRIBUTION'
WRITE (STAT%FHNDL,*) 'PLEASE CHECK THE AMOUNG OF DIRECTIONAL BINS AND THE PRESCRIBED DIRECTIONAL SPREADING'
END IF
ENDDO
END IF
ETOT = 0.
DO ID = 1, NUMDIR
DO IS = 2, NUMSIG
DS = SPSIG(IS) - SPSIG(IS-1)
EAD = 0.5*(SPSIG(IS)*WALOC(IS,ID)+SPSIG(IS-1)*WALOC(IS-1,ID))*DS*DDIR
ETOT = ETOT + EAD
END DO
END DO
IF (PrintLOG) THEN
WRITE (STAT%FHNDL,*) 'HS - INPUTSPECTRA - AFTER 2D', 4.0*SQRT(ETOT)
END IF
ETOT = 0.
EFTOT = 0.
PTAIL_ARR = TAIL_ARR(1) - 1.
ETAIL = 1. / (PTAIL_ARR * (1. + PTAIL_ARR * (FRINTH-1.)))
PTAIL_ARR = TAIL_ARR(1) - 2.
EFTAIL = 1. / (PTAIL_ARR * (1. + PTAIL_ARR * (FRINTH-1.)))
DO ID=1, NUMDIR
DO IS = 1, NUMSIG
OMEG = SPSIG(IS)
EAD = FRINTF * SIGPOW(IS,2) * WALOC(IS,ID)
ETOT = ETOT + EAD
EFTOT = EFTOT + EAD * OMEG
ENDDO
IF (NUMSIG .GT. 3) THEN
EAD = SIGPOW(NUMSIG,2) * WALOC(NUMSIG,ID)
ETOT = ETOT + ETAIL * EAD
EFTOT = EFTOT + EFTAIL * OMEG * EAD
! WRITE (*,*) ETAIL * EAD, EFTAIL * OMEG * EAD, WALOC(NUMSIG,ID)
ENDIF
ENDDO
IF (EFTOT.GT.THR) THEN
TM1 = PI2 * ETOT / EFTOT
ELSE
TM1 = 0.
ENDIF
ETOT = 0.
EFTOT = 0.
PTAIL_ARR = TAIL_ARR(1) - 1.
ETAIL = 1. / (PTAIL_ARR * (1. + PTAIL_ARR * (FRINTH-1.)))
PTAIL_ARR = TAIL_ARR(1) - 3.
EFTAIL = 1. / (PTAIL_ARR * (1. + PTAIL_ARR * (FRINTH-1.)))
DO ID=1, NUMDIR
DO IS=1,NUMSIG
EAD = SIGPOW(IS,2) * WALOC(IS,ID) * FRINTF
OMEG2 = SIGPOW(IS,2)
ETOT = ETOT + EAD
EFTOT = EFTOT + EAD * OMEG2
ENDDO
IF (NUMSIG .GT. 3) THEN
EAD = SIGPOW(NUMSIG,2) * WALOC(NUMSIG,ID)
ETOT = ETOT + ETAIL * EAD
EFTOT = EFTOT + EFTAIL * OMEG2 * EAD
ENDIF
ENDDO
IF (EFTOT .GT. THR) THEN
TM2 = PI2 * SQRT(ETOT/EFTOT)
ELSE
TM2 = 0.
END IF
ETOT = 0.
DO ID = 1, NUMDIR
DO IS = 2, NUMSIG
DS = SPSIG(IS) - SPSIG(IS-1)
EAD = 0.5*(SPSIG(IS)*WALOC(IS,ID)+SPSIG(IS-1)*WALOC(IS-1,ID))*DS*DDIR
ETOT = ETOT + EAD
END DO
END DO
IF (PrintLOG) THEN
WRITE (STAT%FHNDL,*) 'HS - INPUTSPECTRA - AFTER 2D', 4.0*SQRT(ETOT)
WRITE (STAT%FHNDL,*) 'TM01, TM02 & HS', TM1, TM2, 4.0*SQRT(ETOT)
END IF
!FLUSH(DBG%FHNDL)
!FLUSH(STAT%FHNDL)
IF (.FALSE.) THEN ! Write WW3 spectra of the input boundary condition ...
DTHD=360./NUMDIR
RTH0=SPDIR(1)/DDIR
DO ID = 1, NUMDIR
THD(ID)=DTHD*(RTH0+MyREAL(ID-1))
END DO
WRITE (4001,1944) 'WAVEWATCH III SPECTRA', NUMSIG, NUMDIR, 1, 'LAI ET AL'
WRITE (4001,1945) (SPSIG(IS)*INVPI2,IS=1,NUMSIG)
WRITE (4001,1946) (MOD(2.5*PI-SPDIR(ID),PI2),ID=1,NUMDIR)
WRITE (4001,901) 'LAI SPEC', 0., 0., 0., 0., 0., 0., 0.
WRITE (4001,902) ((WALOC(IS,ID)*SPSIG(IS)/PI2*RHOW*G9,IS=1,NUMSIG),ID=1,NUMDIR)
END IF
901 FORMAT ('''',A10,'''',2F7.2,F10.1,2(F7.2,F6.1))
902 FORMAT (7E11.3)
1943 FORMAT ( ' File name : ',A,' (',A,')')
1944 FORMAT ('''',A,'''',1X,3I6,1X,'''',A,'''')
1945 FORMAT (8E10.3)
1946 FORMAT (7E11.3)
END SUBROUTINE
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE SET_WAVE_BOUNDARY
USE DATAPOOL
IMPLICIT NONE
INTEGER :: IP, eIdx
IF (LBCWA .OR. LBCSP) THEN ! Spectrum or parametric boundary condition
IF (LINHOM) THEN
DO IP = 1, IWBMNP
eIdx = IWBNDLC(IP)
AC2(:,:,eIdx) = WBAC(:,:,IP)
END DO
ELSE
DO IP = 1, IWBMNP
eIdx = IWBNDLC(IP)
AC2(:,:,eIdx) = WBAC(:,:,1)
END DO
END IF
END IF
END SUBROUTINE
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE COMPUTE_IFILE_IT(IFILE, IT)
USE DATAPOOL
IMPLICIT NONE
INTEGER, INTENT(OUT) :: IFILE, IT
REAL(rkind) :: DTMP
INTEGER ITMP
DTMP = (MAIN%TMJD-BND_TIME_ALL_FILES(1,1)) * DAY2SEC
ITMP = 0
WRITE(DBG%FHNDL,*)'NUM_NETCDF_FILES_BND =,',NUM_NETCDF_FILES_BND
DO IFILE = 1, NUM_NETCDF_FILES_BND
ITMP = ITMP + NDT_BND_FILE(IFILE)
IF (ITMP .GT. INT(DTMP/SEBO%DELT)) EXIT
END DO
WRITE(DBG%FHNDL,*)'IFILE =,',IFILE
ITMP = SUM(NDT_BND_FILE(1:IFILE-1))
IT = NINT(DTMP/SEBO%DELT) - ITMP + 1
IF (IT .GT. NDT_BND_FILE(IFILE)) THEN
IFILE = IFILE + 1
IT = 1
ENDIF
END SUBROUTINE
!**********************************************************************
!* T. Guerin: Equivalent of COMPUTE_IFILE_IT for WW3 binary files *
!**********************************************************************
SUBROUTINE COMPUTE_IT(IT)
USE DATAPOOL
IMPLICIT NONE
INTEGER, INTENT(OUT) :: IT
REAL(rkind) :: DTMP
INTEGER ITMP
DTMP = (MAIN%TMJD-BND_TIME_ALL_FILES(1,1)) * DAY2SEC
IT = NINT(DTMP/SEBO%DELT) + 1
IF (LBINTER) IT = IT + 1
END SUBROUTINE
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE WAVE_BOUNDARY_CONDITION(WBACOUT)
USE DATAPOOL
IMPLICIT NONE
REAL(rkind), INTENT(OUT) :: WBACOUT(NUMSIG,NUMDIR,IWBMNP)
INTEGER :: IP
!AR: WAVE BOUNDARY
! SPPARM(1): Hs, sign. wave height
! SPPARM(2): Wave period given by user (either peak or mean)
! SPPARM(3): average direction
! SPPARM(4): directional spread
! SPPARM(5): spectral shape (1-4),
! 1 - Pierson-Moskowitz
! 2 - JONSWAP
! 3 - BIN
! 4 - Gauss
! positive peak (+) or mean frequency (-)
! SPPARM(6): directional spreading in degree (1) or exponent (2)
! SPPARM(7): gaussian width for the gauss spectrum 0.1
! SPPARM(8): peak enhancement factor for the JONSWAP spectra 3.
!
! Count number of active boundary points ...
!
IF(LWW3GLOBALOUT) THEN
IF (.NOT. ALLOCATED(WW3GLOBAL)) THEN
ALLOCATE(WW3GLOBAL(8,MNP), stat=istat)
IF (istat/=0) CALL WWM_ABORT('wwm_bdcons, allocate error 7')
END IF
END IF
IF (LBCWA) THEN ! Parametric Wave Boundary is prescribed
IF (PrintLOG) THEN
WRITE(STAT%FHNDL,'("+TRACE...",A)') 'Parametric Wave Boundary Condition is prescribed'
END IF
IF (LINHOM) THEN ! Inhomogenous in space
IF (LBCSE) THEN ! Unsteady in time
SPPARM = 0.
IF (IBOUNDFORMAT == 1) THEN ! WWM
CALL READWAVEPARWWM
ELSE IF (IBOUNDFORMAT == 2) THEN ! WW3
#ifdef NCDF
CALL READ_NETCDF_WW3_PARAM
#else
CALL WWM_ABORT('compile with netcdf for IBOUNDFORMAT=2')
#endif
CALL INTER_STRUCT_BOUNDARY(NDX_BND,NDY_BND,DX_BND,DY_BND,OFFSET_X_BND,OFFSET_Y_BND,SPPARM)
IF (LWW3GLOBALOUT) CALL INTER_STRUCT_DOMAIN(NDX_BND,NDY_BND,DX_BND,DY_BND,OFFSET_X_BND,OFFSET_Y_BND,WW3GLOBAL)
ELSE IF (IBOUNDFORMAT == 4) THEN ! WWM SPPARM netcdf file
#ifdef NCDF
CALL READ_NETCDF_BOUNDARY_SPPARM
#else
CALL WWM_ABORT('compile with netcdf for IBOUNDFORMAT=4')
#endif
ELSE IF (IBOUNDFORMAT == 5) THEN ! WAM format of waves
CALL WWM_ABORT('No possibility of using parametric boundary for IBOUNDFORMAT=5 yet')
ELSE
CALL WWM_ABORT('Wrong value of IBOUNDFORMAT')
END IF
ELSE ! Steady ...
SPPARM = 0.
IF (IBOUNDFORMAT == 1) THEN
CALL READWAVEPARWWM
ELSE IF (IBOUNDFORMAT == 2) THEN
#ifdef NCDF
CALL READ_NETCDF_WW3_PARAM
#else
CALL WWM_ABORT('compile with netcdf for IBOUNDFORMAT=2')
#endif
CALL INTER_STRUCT_BOUNDARY(NDX_BND,NDY_BND,DX_BND,DY_BND,OFFSET_X_BND,OFFSET_Y_BND,SPPARM)
END IF
END IF ! LBCSE ...
DO IP = 1, IWBMNP
CALL SPECTRAL_SHAPE(SPPARM(:,IP),WBACOUT(:,:,IP),.FALSE.,'CALL FROM WB 1', USE_OPTI_SPEC_SHAPE_BOUC)
END DO
ELSE ! Homogenous in space
IF (IWBMNP .gt. 0) THEN
IF (LBCSE) THEN ! Unsteady in time
IF (IBOUNDFORMAT == 1) THEN
CALL READWAVEPARWWM
ELSE IF (IBOUNDFORMAT == 2) THEN
#ifdef NCDF
CALL READ_NETCDF_WW3_PARAM
#else
CALL WWM_ABORT('compile with netcdf for IBOUNDFORMAT=2')
#endif
CALL INTER_STRUCT_BOUNDARY(NDX_BND,NDY_BND,DX_BND,DY_BND,OFFSET_X_BND,OFFSET_Y_BND,SPPARM)
IF (LWW3GLOBALOUT) CALL INTER_STRUCT_DOMAIN(NDX_BND,NDY_BND,DX_BND,DY_BND,OFFSET_X_BND,OFFSET_Y_BND,WW3GLOBAL)
CALL SPECTRAL_SHAPE(SPPARM(:,1),WBACOUT(:,:,1), .FALSE.,'CALL FROM WB 3', USE_OPTI_SPEC_SHAPE_BOUC)
ELSE IF (IBOUNDFORMAT == 4) THEN
#ifdef NCDF
CALL READ_NETCDF_BOUNDARY_SPPARM
#else
CALL WWM_ABORT('compile with netcdf for IBOUNDFORMAT=4')
#endif
ELSE IF (IBOUNDFORMAT == 5) THEN
CALL WWM_ABORT('No possibility of using parametric boundary for IBOUNDFORMAT=5')
ELSE
CALL WWM_ABORT('No possibility of using parametric boundary condition')
END IF
ELSE ! Steady in time ...
SPPARM = 0.
IF (LMONO_IN) THEN
SPPARM(1,1) = WBHS * SQRT(2.)
ELSE
SPPARM(1,1) = WBHS
END IF
SPPARM(2,1) = WBTP
SPPARM(3,1) = WBDM
SPPARM(4,1) = WBDS
SPPARM(5,1) = WBSS
SPPARM(6,1) = WBDSMS
SPPARM(7,1) = WBGAUSS
SPPARM(8,1) = WBPKEN
CALL SPECTRAL_SHAPE(SPPARM(:,1),WBACOUT(:,:,1),.FALSE.,'CALL FROM WB 4', USE_OPTI_SPEC_SHAPE_BOUC)
END IF ! LBCSE
END IF
END IF ! LINHOM
END IF
IF (LBCSP) THEN ! Spectrum is prescribed
IF (LINHOM) THEN ! The boundary conditions is not homogenous!
IF (LBSP1D) CALL WWM_ABORT('No inhomogenous 1d spectra boundary cond. available')
IF (IBOUNDFORMAT == 1) THEN ! WWM
!CALL READSPEC2D
CALL WWM_ABORT('No inhomogenous 2d spectra boundary cond. available in WWM Format')
END IF
IF (IBOUNDFORMAT == 2) THEN ! WW3 KM
CALL GET_BINARY_WW3_SPECTRA(WBACOUT)
END IF
IF (IBOUNDFORMAT == 4) THEN ! WWM WBAC netcdf
#ifdef NCDF
CALL READ_NETCDF_BOUNDARY_WBAC(WBACOUT)
#else
CALL WWM_ABORT('compile with netcdf for IBOUNDFORMAT=4')
#endif
END IF
IF (IBOUNDFORMAT == 5) THEN ! WWM WBAC netcdf
#ifdef GRIB_API_ECMWF
CALL READ_GRIB_WAM_BOUNDARY_WBAC(WBACOUT)
#else
CALL WWM_ABORT('compile with GRIB_API for IBOUNDFORMAT=5')
#endif
END IF
ELSE ! The boundary conditions is homogeneous in space !
IF (LBSP1D) THEN ! 1-D Spectra is prescribed
IF (PrintLOG) THEN
WRITE(STAT%FHNDL,'("+TRACE...",A)') '1d Spectra is given as Wave Boundary Condition'
END IF
CALL READSPEC1D(LFIRSTREAD)
CALL SPECTRUM_INT(WBACOUT(:,:,1))
ELSE IF (LBSP2D) THEN ! 2-D Spectra is prescribed
IF (PrintLOG) THEN
WRITE(STAT%FHNDL,'("+TRACE...",A)') '2d Spectra is given as Wave Boundary Condition'
END IF
IF (IBOUNDFORMAT == 1) THEN
CALL READSPEC2D(LFIRSTREAD)
ELSE IF (IBOUNDFORMAT == 2) THEN !KM
CALL GET_BINARY_WW3_SPECTRA(WBACOUT)
END IF
!AR: check what is goin on there ...
!CALL SPECTRUM_INT(WBACOUT(:,:,1)) !KM
END IF ! LBSP1D .OR. LBSP2D
END IF ! LINHOM
ENDIF ! LBCSP
END SUBROUTINE
!**********************************************************************
!* *
!**********************************************************************
SUBROUTINE SET_WAVE_BOUNDARY_CONDITION
USE DATAPOOL
IMPLICIT NONE
CHARACTER(len=29) :: CHR
IF (LNANINFCHK) THEN
WRITE(DBG%FHNDL,*) ' ENTERING SET_WAVE_BOUNDARY_CONDITION ', SUM(AC2)
IF (SUM(AC2) .NE. SUM(AC2)) CALL WWM_ABORT('NAN IN BOUNDARY CONDTITION l.1959')
ENDIF
IF ( MAIN%TMJD > SEBO%TMJD-1.E-8 .AND. MAIN%TMJD < SEBO%EMJD ) THEN ! Read next time step from boundary file ...