-
Notifications
You must be signed in to change notification settings - Fork 0
/
Copy pathCSP_PHENOL.for
1237 lines (1009 loc) · 45.8 KB
/
CSP_PHENOL.for
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
!======================================================================
!
! CSP_PHENOL Subroutine, based on PHENOL_SC, Subroutine, O.H. Daza
!
! PURPOSE: Simualte phenological development of the sugarcane crop
!----------------------------------------------------------------------
! REVISION HISTORY
! 03/10/01 OHD wrote it
! 07/27/03 FSR reformatted to FORTRAN 77
! 08/22/03 FSR integrated into DSSAT 4.0
! 07/26/2004 CHP Removed variables which were not being used
! 07/07/2005 FSR Modified
!----------------------------------------------------------------------
! Called from: CASUPRO
! Calls: CSP_INPHENOL in CSP_INPHENOL.FOR
! CSP_IPPHENOL in CSP_IPPHENOL.FOR
! CURV in UTILS.FOR
!----------------------------------------------------------------------
!
SUBROUTINE CSP_PHENOL(CONTROL, ISWITCH, FILECC,
& DTPI, Kill, LI, M, SOILPROP, TGROAV, !Input
& TURFAC, XLFNUM, XStkNum, YLfFac, YLFSZ, YRDOY, !Input
& YRPLT, YRSIM, !Input
& CropTypeCode, DeltaLeafArea, DeltaLeafNum, !Output
& DLFN, DPLARF, DRPP, DTX, TillerCount, !Output
& LeafNum, MinGr, MDATE, !Output
& NVEG0, PhenoStage, PHTHRS, Ph1P, PI1, PI2, ROWSPC, !Output
& Smax, StalkState, STNAME, STGDOY, StkHrNO, !Output
& SumTTD, SumTTG, SumTTStalk, TDUMX, VSTAGE, !Output
& XLAI, YREMRG) !Output
!----------------------------------------------------------------------
USE ModuleDefs !Definitions of constructed variable types,
! which contain control information, soil
! parameters, hourly weather data.
IMPLICIT NONE
SAVE
!----------------------------------------------------------------------
! Variable specification
CHARACTER*1 ISIMI, ISWWAT, IDETG, PLME
CHARACTER*2 CROP
!CHARACTER*3 CTMP(4), DLTYP(4)
CHARACTER*4 StalkState(NumOfStalks,10)
CHARACTER*10 STNAME(20)
CHARACTER(11) :: CropTypeName
CHARACTER*13 OUTPHNL
CHARACTER*30 FILEIO
CHARACTER*92 FILECC
INTEGER Cond, EndStalk, ERRNUM, FROP, DYNAMIC, DOY, DAS
INTEGER DAP, TillerCount, MDATE, NLAYR, NOUTPHN, NVEG0, RUN
INTEGER Smax, TIMDIF, YEAR, YREMRG, YRPLT, YRSIM, YRDOY
INTEGER livecount, Day, Phase, Stalk
INTEGER CropTypeCode
INTEGER PhenoStage, NewStalk
INTEGER Temp , WLUN
INTEGER, PARAMETER :: NumOfTemp = 5, NumOfLeaves = 40
INTEGER, PARAMETER :: NumOfStages = 5, NumOfPhases = 4
INTEGER, DIMENSION(NumOfStages) :: DayOfStage,
& DaysAfterPlantOfStage, STGDOY
INTEGER, DIMENSION(1:NumOfStalks) :: Kill
REAL LI, LI1, VT, RTR, XLI(7), YVTR(7)
REAL XLFNUM(7), YLFSZ(7), XStkNum(9), YLfFac(9)
REAL TABEX ! Function subroutine - Lookup utility
REAL TURFAC, DRPP, DTX
REAL SDEPTH
REAL PLANTS, PLTPOP, ROWSPC
REAL TDUMX, ExcsTTD
REAL CURV ! Function subroutine
REAL Tbase, TbaseStalk
REAL M, MinGr
REAL Ph1P, Ph1R, Ph2, Ph3, Ph4
REAL DTPI, PI1, PI2
REAL So, Go, Gmax
REAL DepthRateOfEmer, DeltaDepthOfEmer, DepthToEmer
REAL TGROAV, DeltaTTD, SumTTD, DeltaTTG
REAL RTNFAC, StkHrNO, SumTTG, TELOM
REAL XLAI
REAL VSTAGE, EnviroFactor
REAL GrowFrac, GrowTime, RipeFrac, RipeTime
LOGICAL FEXIST !, FIRST
REAL :: DPLARF, DLFN, NewTiller, TillerExcess
REAL, DIMENSION(NL) :: LL, DUL, SAT, DLAYR ! , SW, ST
REAL, DIMENSION(0:NumOfPhases) :: FNSTR, FPSTR, FSW, FT, FUDAY,
& WSENP, NSENP
REAL, DIMENSION(1:NumOfStalks) :: DeltaTillerNum
REAL, DIMENSION(0:NumOfDays,NumOfStalks) :: DeltaLeafNum, LeafNum,
& DeltaLeafArea
REAL, DIMENSION(0:NumOfDays,NumOfStalks) :: SumTTStalk
REAL, DIMENSION(NumOfPhases) :: PHTHRS, PHZACC
REAL, DIMENSION(NumOfPhases) :: OptStageDur, CumOptStageDur
REAL, DIMENSION(NumOfTemp) :: TB, TO1, TO2, TM
!-----------------------------------------------------------------------
! Define constructed variable types based on definitions in
! ModuleDefs.for.
! The variable "CONTROL" is of type "ControlType".
TYPE (ControlType) CONTROL
! The variable "SOILPROP" is of type "SoilType".
TYPE (SoilType) SOILPROP
! The variable "ISWITCH" is of type "SwitchType".
TYPE (SwitchType) ISWITCH
! No output for fallow crop (or if IDETG = N? FSR)
CROP = CONTROL % CROP
IDETG = ISWITCH % IDETG
! Transfer values from constructed data types into local variables.
DAS = CONTROL % DAS
DYNAMIC = CONTROL % DYNAMIC
FILEIO = CONTROL % FILEIO
FROP = CONTROL % FROP
RUN = CONTROL % RUN
YRDOY = CONTROL % YRDOY
YRSIM = CONTROL % YRSIM
DLAYR = SOILPROP % DLAYR
DUL = SOILPROP % DUL
LL = SOILPROP % LL
NLAYR = SOILPROP % NLAYR
SAT = SOILPROP % SAT
ISWWAT = ISWITCH % ISWWAT
!-----------------------------------------------------------------------
! Start execution part
!-----------------------------------------------------------------------
! Days after start of simulation
! DAS = MAX(0,TIMDIF(YRSIM,YRDOY))
! Days after planting
DAP = MAX(0,TIMDIF(YRPLT,YRDOY))
! Bud planting density is read in FILEX as PLANTS; may want to use this
! method in the future, since it corresponds to production practice
! BudPlantDensity = ((1 / BudSpacing) + 1) / RowSpacing
! Right
! BudPlantDensity = 1 / (BudSpacing * RowSpacing)
!**********************************************************************
!**********************************************************************
! Run Initialization - Called once per simulation
!**********************************************************************
IF (DYNAMIC .EQ. RUNINIT) THEN
DPLARF = 0 !Leaf area for the whole plant
EndStalk = 0
NewTiller = 0.0
RTR = 0.0
LI = 0.0 ! Initialize here since Phenol is called before Photo (FSR)
TillerCount = 0.0
!----------------------------------------------------------------------
!FSR - Introduced another output file that is readable by G-Build.
! Included the following line in UTILS.FOR:
! CASE ('OUTPHNL'); LUN = 48 !PHENOLOGY.OUT - FSR
! IF (IDETG .EQ. 'Y') THEN !Whenever PlantGrow.Out is selected
IF (INDEX('Y' ,ISWITCH%IDETG) > 0 .AND.
& INDEX('YDA',ISWITCH%IDETL) > 0) THEN
OUTPHNL = 'CSP_PHEN.OUT'
CALL GETLUN('OUTPHNL', NOUTPHN)
ENDIF
!----------------------------------------------------------------------
!----------------------------------------------------------------------
! Subroutine CSP_IPPHENOL reads required phenology variables from
! input files. Subroutine INPHENOL initializes data for beginning
! of run (Once only initialization).
!----------------------------------------------------------------------
CALL CSP_IPPHENOL(CONTROL, FILECC,
& CROP, DTPI, Gmax, Go, ISIMI, !Output
& LI1, MinGr, Ph1P, Ph1R, Ph2, Ph3, !Output
& Ph4, PI1, PI2, PLANTS, PLME, PLTPOP, !Output
& RTNFAC, ROWSPC, SDEPTH, Smax, So, !Output
& StkHrNO, TB, TELOM, TM, TO1, !Output
& TO2, XLFNUM, XLI, XStkNum, YLfFac, !Output
& YLFSZ, YVTR) !Output
! PHTHRS calculate this here and make it input in CSP_INPHENOL
! It will become 2 or greater after the first run of simulation
! see MESIC and others
! CropTypeCode Crop type code (1 : plant crop; >=2 : ratoon crop)
IF (PLME .EQ. 'R') THEN !PLME Planting Method, is set in FILEX.
CropTypeCode = 2 ! 1 indicates plant cane, and 2 is the
! first ratoon option. Second and higher
! ratoons will require further work
ELSE
CropTypeCode = 1
END IF
IF (CropTypeCode == 1) THEN
PHTHRS(1) = Ph1P
CropTypeName = "Plant crop"
! END IF
ELSE IF (CropTypeCode >= 2) THEN
PHTHRS(1) = Ph1R
CropTypeName = "Ratoon crop"
END IF
! PHTHRS from "phase threshold" (FSR)
PHTHRS(2) = SDEPTH * 10 / Ph2 ! Determine PhenoStages
PHTHRS(3) = Ph3 ! Ph1, 2, 3 & 4 parameters in .ECO
PHTHRS(4) = Ph4
CALL GETLUN('WORK.OUT', WLUN)
OPEN(UNIT = WLUN, FILE = "WORK.OUT", STATUS = "UNKNOWN",
& ACTION = "WRITE", POSITION = "APPEND")
WRITE(WLUN,'(1X,"RESULTS FROM CSP_IPPHENOL.FOR")')
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"Input file: ",A)') FILEIO
WRITE(WLUN,'(1X,"Crop : ",A2)') CROP
WRITE(WLUN,'(1X,"ISIMI : ",A1)') ISIMI
WRITE(WLUN,'(1X,"ISWWAT: ",A1)') ISWWAT
WRITE(WLUN,*)
WRITE(WLUN,'(1x,"Plant population at seeding:
& ",F8.1," plants/m²")') PLANTS
WRITE(WLUN,'(1x,"Plant population at emergence:
& ",F8.1," plants/m²")') PLTPOP
WRITE(WLUN,'(1X,"PLME : ",A1)') PLME
WRITE(WLUN,'(1x,"Row spacing :
& ",F8.1," cm")') ROWSPC
WRITE(WLUN,'(1x,"Planting depth :
& ",F8.1," cm")') SDEPTH
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"Output from cultivar file .CUL")')
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"Phyllochron interval 1:",F8.1," °C-day")') PI1
WRITE(WLUN,'(1X,"Phyllochron interval 2:",F8.1," °C-day")') PI2
WRITE(WLUN,'(1X,"Phyllochron interval :",F8.1," °C-day")') DTPI
WRITE(WLUN,*)
!!!WRITE(WLUN,'(1X,"Smax:",F8.1," # stalks/stubble")') Smax
WRITE(WLUN,'(1X,"Smax:",I8," # stalks/stubble")') Smax
WRITE(WLUN,'(1X,"So :",F8.1," # stalks/stubble")') So
WRITE(WLUN,'(1X,"Gmax:",F8.1," °C-day")') Gmax
WRITE(WLUN,'(1X,"Go :",F8.1," °C-day")') Go
WRITE(WLUN,*)
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"Phase 1 - plant crop :", F8.1, " ºC-day")') Ph1P
WRITE(WLUN,'(1X,"Phase 1 - ratoon crop:", F8.1, " ºC-day")') Ph1R
WRITE(WLUN,'(1X,"Phase 2 :
& ", F8.1, " mm/(ºC-day)")') Ph2
WRITE(WLUN,'(1X,"Phase 3 :", F8.1, " ºC-day")') Ph3
WRITE(WLUN,'(1X,"Phase 4 :", F8.1, " ºC-day")') Ph4
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"Phase PHTHRS ºC-day")')
DO Phase = 1, NumOfPhases
WRITE(WLUN,'(1X,I5,1X,F7.1)') Phase, PHTHRS(Phase)
END DO
WRITE(WLUN,*)
WRITE(WLUN,'(1X," TB TO1 TO2 TM ºC")')
DO Temp = 1,NumOfTemp
WRITE(WLUN,'(1X, 4F5.1)') TB(Temp), TO1(Temp), TO2(Temp),
& TM(Temp)
END DO
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"Phase WSENP NSENP")')
DO Phase = 1, NumOfPhases
WRITE(WLUN,'(1X,I5,2(1X,F5.1))') Phase, WSENP(Phase),
& NSENP(Phase)
END DO
! WRITE(WLUN,'(1X,"Crop: ", A)') CropTypeName(:6)
!----------------------------------------------------------------------
! Subroutine CSP_INPHENOL calculates required phenology variables
! from values read in input files. Subroutine CSP_INPHENOL initializes
! data for beginning of run (Once only initialization).
!----------------------------------------------------------------------
CALL CSP_INPHENOL(
& CROP, PHTHRS, TB, TO1, !Input
& CumOptStageDur, OptStageDur, STNAME) !Output
! TEMPORARY: Statements to test output from CSP_INPHENOL above
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"RESULTS FROM CSP_INPHENOL.FOR")')
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"Phase OptStageDur CumOptStageDur STNAME")')
WRITE(WLUN,'(1X," (days) (days)")')
DO Phase = 1, NumOfPhases
WRITE(WLUN,'(1X,I5,1X,F11.1,1X,F14.1,2X,A)') Phase,
& OptStageDur(Phase), CumOptStageDur(Phase),
& STNAME(Phase)
END DO
! Base temperatures
Tbase = TB(1)
TbaseStalk = TB(2)
! This is an input value from file .CUL (read by CSP_IPPHENOL)
DepthRateOfEmer = Ph2 ! e.g., 0.8 mm (soil depth) / (ºC-day)
! For P module (FSR):
GrowFrac = 0.0
RipeFrac = 0.0
GrowTime = PHTHRS(1) + PHTHRS(2) + PHTHRS(3)
RipeTime = PHTHRS(4)
!**********************************************************************
!**********************************************************************
! Seasonal initialization - run once per season
!**********************************************************************
ELSE IF (DYNAMIC .EQ. SEASINIT) THEN
!-----------------------------------------------------------------------
IF (INDEX('Y' ,ISWITCH%IDETG) > 0 .AND.
& INDEX('YDA',ISWITCH%IDETL) > 0) THEN
!FSR - Initialize daily phenology.out file
INQUIRE (FILE = OUTPHNL, EXIST = FEXIST)
IF (FEXIST) THEN
OPEN (UNIT = NOUTPHN, FILE = OUTPHNL, STATUS = 'OLD',
& IOSTAT = ERRNUM, POSITION = 'APPEND') ! Was 'APPEND' FSR
!!! FIRST = .FALSE.
ELSE
OPEN (UNIT = NOUTPHN, FILE = OUTPHNL, STATUS = 'NEW',
& IOSTAT = ERRNUM)
WRITE(NOUTPHN,'("*PHENOLOGY ASPECTS OUTPUT FILE")')
!!! FIRST = .TRUE.
ENDIF
!Write headers
! - - - - - New CASUPRO header v v v v - - - - - - -
CALL HEADER(SEASINIT, NOUTPHN, RUN)
WRITE(NOUTPHN,'("@YEAR DOY DAS DAP PHSTG NSKST STTD",
& " DTTG STTG XLAI")', ADVANCE="NO")
! Header line is complicated to keep column headings same length
! by inserting underscore where stalk number < 10. FSR
IF (Smax <= 9) THEN
DO Stalk = 1,Smax
WRITE(NOUTPHN,'(3X"LN_",I1)',ADVANCE="NO") Stalk
END DO
ELSE
DO Stalk = 1,9
WRITE(NOUTPHN,'(3X"LN_",I1)',ADVANCE="NO") Stalk
END DO
DO Stalk = 10,Smax
WRITE(NOUTPHN,'(3X"LN",I2)',ADVANCE="NO") Stalk
END DO
ENDIF
IF (Smax <= 9) THEN
DO Stalk = 1,Smax
WRITE(NOUTPHN,'(2X"DLAR_",I1)',ADVANCE="NO") Stalk
END DO
ELSE
DO Stalk = 1,9
WRITE(NOUTPHN,'(2X"DLAR_",I1)',ADVANCE="NO") Stalk
END DO
DO Stalk = 10,Smax
WRITE(NOUTPHN,'(2X"DLAR",I2)',ADVANCE="NO") Stalk
END DO
ENDIF
ENDIF
!-----------------------------------------------------------------------
! Ph1P Threshold to sprouting - Plant cane, °C-day
! Ph1R Threshold to sprouting - Ratoon cane, °C-day
!-----------------------------------------------------------------------
IF (CropTypeCode >= 2) THEN
PHTHRS(1) = Ph1R
CropTypeName = "Ratoon crop"
! YRPLT = YRRAT ! Date of ratooning
END IF
PhenoStage = 0 ! Added by CHP and FSR to solve initialization
VSTAGE = 0.0 ! problems, leading to emergence a day earlier
DeltaTTD = 0.0 ! in batch runs. 10-08-2003
DeltaTTG = 0.0
ExcsTTD = 0
Stalk = 0
NewStalk = 0
SumTTD = 0 ! Summation of thermal time for leaves
SumTTG = 0 ! Summation of thermal time for stalks
DepthToEmer = 0
NVEG0 = 10000
! These can be incorporated later on to affect development. OHD
DRPP = 0.0 ! Photoperiod days which occur in a real day
! (photoperiod days / day)
DTX = 0.0 ! Thermal time that occurs in a real day based on
! vegetative development temperature function
! (thermal days / day)
TDUMX = 0.0 ! Photo-thermal time that occurs in a real day
! based on early reproductive development
! temperature function
XLAI = 0.0 ! Leaf Area Index
! Initialization of arrays for each phase
DO Phase = 1, NumOfPhases
FNSTR(Phase) = 1. ! stress-by-phase not yet implemented
FPSTR(Phase) = 1. ! check definitions at end of file - FSR
FSW(Phase) = 1.
FT(Phase) = 0.
FUDAY(Phase) = 0.
PHZACC(Phase) = 0. ! Of this list, only this is currently used
END DO
DO Stalk = 1, Smax
DeltaTillerNum(Stalk) = 0.
END DO
! Cond 1 - LIVE or DEAD; 2 - PRIM or TILR; 3 - 10 not yet used
DO Stalk = 1, Smax
DO Cond = 1, 10
StalkState(Stalk,Cond) = ' '
END DO
END DO
DO Day = 0, NumOfDays
! DO Day = 1, NumOfDays
! Initialization of variables for stalk appearance and leaf area
DO Stalk = 1, Smax
LeafNum(Day,Stalk) = 0.
SumTTStalk(Day,Stalk) = 0.!Thermal time for stalk appearance
DeltaLeafArea(Day,Stalk)= 0.
DeltaLeafNum(Day,Stalk) = 0.
END DO ! Stalk
END DO ! Day
! TEMPORARY: Statements to test output from PHENOL_SC
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"RESULTS FROM CSP_PHENOL.FOR - SEASINIT")')
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"YRSIM:",I7)') YRSIM
WRITE(WLUN,'(1X,"YRPLT:",I7)') YRPLT
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"CropTypeName: ",A11)') CropTypeName
WRITE(WLUN,*)
WRITE(WLUN,'(1X,"RESULTS FROM CSP_PHENOL.FOR - INTEGR")')
WRITE(WLUN,*)
WRITE(WLUN,'(1X," YRDOY DAS DAP PHSTG NSKST STTD
&DTTG STTG XLAI")', ADVANCE="NO")
! DO Stalk = 1, Smax
! WRITE(WLUN,'(" LN(",I2,")")',ADVANCE="NO") Stalk
! END DO
DO Stalk = 1,Smax
WRITE(WLUN,'(3X"LN",I2)',ADVANCE="NO") Stalk
END DO
DO Stalk = 1, Smax
WRITE(WLUN,'(2X"DLAR",I2)',ADVANCE="NO") Stalk
END DO
!***********************************************************************
!***********************************************************************
! Daily Rate calculations
!***********************************************************************
ELSE IF (DYNAMIC .EQ. RATE) THEN
!-----------------------------------------------------------------------
! Compute temperature, daylength, and water effects on development.
! Daylength and water effects on development not included yet.
!-----------------------------------------------------------------------
!************************************************************************
! CODE per Fernando Villegas FVT / CHP 9/16/2004
!************************************************************************
IF (YRDOY > YRPLT) THEN
! Calculate thermal time in the day for primary stalk & leaf appearance
IF (TGROAV < TB(1) .OR. TGROAV > TM(1)) THEN
DeltaTTD = 0.
ELSEIF (TGROAV >= TB(1) .AND. TGROAV < TO1(1)) THEN
DeltaTTD = TGROAV - TB(1)
ELSEIF (TGROAV >= TO1(1) .AND. TGROAV <= TO2(1)) THEN
DeltaTTD = TO1(1) - TB(1)
ELSEIF (TGROAV > TO2(1) .AND. TGROAV <= TM(1)) THEN
DeltaTTD = (TO1(1) - TB(1))*(TM(1) - TGROAV)/(TM(1) - TO2(1))
ENDIF
! Calculate thermal time in the day for tiller appearance
IF (TGROAV < TB(2) .OR. TGROAV > TM(2)) THEN
DeltaTTG = 0.
ELSEIF (TGROAV >= TB(2) .AND. TGROAV < TO1(2)) THEN
DeltaTTG = TGROAV - TB(2)
ELSEIF (TGROAV >= TO1(2) .AND. TGROAV <= TO2(2)) THEN
DeltaTTG = TO1(2) - TB(2)
ELSEIF (TGROAV > TO2(2) .AND. TGROAV <= TM(2)) THEN
DeltaTTG = (TO1(2) - TB(2))*(TM(2) - TGROAV)/(TM(2) - TO2(2))
ENDIF
ENDIF
! DTX Thermal time that occurs in a real day based on vegetative
! development temperature function (thermal days / day)
! It could be substituted in this module. (FSR)
! TURFAC Water stress factor for expansion (0 - 1)
! DTX could be substituted in
DTX = CURV('lin', Tb(1), To1(1), To2(1), Tm(1), TGROAV)
!************************************************************************
! Environmental factor as a function of water stress factor
! This factor is taken from the equation below defined in PHENOL.FOR and adapted
! for the sugarcane model. Will P-stress be added here? FSR
! VSTAGE = VSTAGE + DTX * TRIFOL * EVMOD * TURFAC * (1.0 - XPOD)
EnviroFactor = TURFAC
! EnviroFactor = TURFAC*DTX !Use DTX = DeltaTTD/(TO1(1)-TB(1))
! Note: DTX is replaced with DeltaTTD.
! TRIFOL not used in CASUPRO, is the rate of appearance of leaves
! on mainstem; replaced with DeltaLeafNum(DAS, Stalk).
! EVMOD Modifies rate of development; not used in CASUPRO.
! XPOD Growth partitioning to pods.
!
! So only TURFAC remains, & use of EnviroFactor needs to be evaluated.
!************************************************************************
!FSR added for Plant P - 08/25/2006
IF (GrowTime > 0) THEN
GrowFrac = (PHZACC(1) + PHZACC(2) + PHZACC(3)) / GrowTime
ELSE
GrowFrac = 0.0
ENDIF
IF (RipeTime > 0) THEN
RipeFrac = PHZACC(4) / RipeTime
ELSE
RipeFrac = 0.0
ENDIF
!************************************************************************
SELECT CASE (PhenoStage) ! Checks PhenoStage value
!=============
! CASE (1) ! planting - sprouting phase
!=============
!=============
CASE (2) ! sprouting - emergence phase
!=============
IF (CropTypeCode == 1) THEN
! Bud planting depth in cm
IF (DeltaTTD > 0) THEN
DeltaDepthOfEmer = DepthRateOfEmer * DeltaTTD / 10 !mm to cm
END IF
END IF
!=============
CASE (3) ! emergence - stalk growth phase
!=============
!******************* Stalk development code (FSR) **********************
IF (LI < LI1) THEN !Low competition among tillers
!! Use TELOM for tiller emergence thermal time (°C-day):.
RTR = 1/TELOM
ELSE ! Calculate Relative Tillering Rate following
! Bezuidenhout et al 2003
RTR = MAX(0.0,TABEX(YVTR, XLI, LI, 6))
END IF ! (LI < LI1)
DO Stalk = 1, Smax
IF (StalkState(Stalk,1) .EQ. 'LIVE') THEN
DeltaTillerNum(Stalk) = DeltaTTG * RTR
ELSE ! Stalk > NewStalk (not alive)
DeltaTillerNum(Stalk) = 0
END IF ! (StalkState(Stalk,1) .EQ. 'LIVE')
END DO
!***********************************************************************
! Section for leaf development
!-----------------------------------------------------------------------
! Calculates the rate of increase of leaf number on each stalk
DO Stalk = 1, Smax
IF (StalkState(Stalk,1) .EQ. 'LIVE') THEN
IF (SumTTStalk(DAS - 1, Stalk) <= DTPI) THEN
DeltaLeafNum(DAS, Stalk) = DeltaTTD / PI1
ELSE
DeltaLeafNum(DAS, Stalk) = DeltaTTD / PI2
END IF
ELSE
DeltaLeafNum(DAS, Stalk) = 0
END IF
!Correction due to temperature and water deficit
DeltaLeafNum(DAS, Stalk) = DeltaLeafNum(DAS, Stalk)
& * EnviroFactor
!--- New look-up code replaces former M1,M2 leaf area parameters ------
M = MAX(TABEX(YLFSZ, XLFNUM, LeafNum(DAS-1, Stalk), 7), 1.E-06)
& * MAX(TABEX(YLfFac, XStkNum, Float(Stalk), 9), 1.E-06)
DeltaLeafArea(DAS,Stalk) = M * DeltaLeafNum(DAS,Stalk)
! NOTE: MAX function avoids dividing by zero in CSP_Grow_Cane
! when LeafNum = 0, M = 0
!-----------------------------------------------------------------------
END DO
!=============
CASE (4) ! stalk growth - flowering phase
!=============
!******************* Stalk development code (FSR) **********************
! Calculate Relative Tillering Rate based on Bezuidenhout et al 2003
RTR = MAX(0.0,TABEX(YVTR, XLI, LI, 6))
DO Stalk = 1, Smax
IF (StalkState(Stalk,1) .EQ. 'LIVE') THEN
DeltaTillerNum(Stalk) = DeltaTTG * RTR
ELSE ! Stalk > NewStalk (not alive)
DeltaTillerNum(Stalk) = 0
END IF ! (StalkState(Stalk,1) .EQ. 'LIVE')
END DO
!***********************************************************************
! Section for leaf development
DO Stalk = 1, Smax
!-----------------------------------------------------------------------
! Calculates the rate of increase of leaf number in each stalk
IF (StalkState(Stalk,1) .EQ. 'LIVE') THEN
IF (SumTTStalk(DAS - 1, Stalk) <= DTPI) THEN
DeltaLeafNum(DAS, Stalk) = DeltaTTD / PI1
ELSE
DeltaLeafNum(DAS, Stalk) = DeltaTTD / PI2
END IF
ELSE
DeltaLeafNum(DAS, Stalk) = 0
END IF
!Correction due to temperature and water deficit
DeltaLeafNum(DAS, Stalk) = DeltaLeafNum(DAS, Stalk)
& * EnviroFactor
!--- New look-up code replaces former M1,M2 leaf area parameters ------
M = TABEX(YLFSZ, XLFNUM, LeafNum(DAS-1, Stalk), 7)
& * TABEX(YLfFac, XStkNum, float(Stalk), 9)
DeltaLeafArea(DAS,Stalk) = M * DeltaLeafNum(DAS,Stalk)
! NOTE: MAX function avoids dividing by zero in CSP_Grow_Cane
! when LeafNum = 0, M = 0
!-----------------------------------------------------------------------
END DO
END SELECT
!----------------------------------------------------------------
! Increment of leaf number in main stalk
DLFN = DeltaLeafNum(DAS,1)
!**********************************************************************
!**********************************************************************
! Daily Integration
!**********************************************************************
ELSE IF (DYNAMIC .EQ. INTEGR) THEN
!**********************************************************************
DTX = CURV('lin', Tb(1), To1(1), To2(1), Tm(1), TGROAV)
! DTX calculation repeated here since RATE section does not exist
! in CSP_ROOTS, and value of DTX was not being passed.
IF (YRDOY .EQ. YRPLT) THEN ! It will be the same for ratooning
PhenoStage = 1 ! Planting or ratooning is set
STGDOY(1) = YRDOY
DayOfStage(1) = YRDOY ! Sets the day
DaysAfterPlantOfStage(1) = DAP ! Sets de number of days
! DeltaTTD = 0 ! Thermal time is not accounted the first day of simulation
END IF
! Integrates thermal time for days greater than day of planting
! For ratooning YRPLT = YRRAT - Make this later!
IF (YRDOY > YRPLT) THEN
!Integrates thermal time only for DeltaTTD > 0
IF (DeltaTTD > 0) SumTTD = SumTTD + DeltaTTD
IF ((PhenoStage == 3) .OR. (PhenoStage == 4)) THEN
!Integrates thermal time only for DeltaTTG > 0
IF (DeltaTTG > 0) SumTTG = SumTTG + DeltaTTG
END IF
END IF
SELECT CASE (PhenoStage) ! Checks PhenoStage value
!=============
CASE (1)
!=============
IF (CropTypeCode == 1) THEN
PHZACC(1) = PHZACC(1) + DeltaTTD
IF (PHZACC(1) > PHTHRS(1)) THEN
PhenoStage = 2
STGDOY(2) = YRDOY
DayOfStage(2) = YRDOY
DaysAfterPlantOfStage(2) = DAP
ExcsTTD = PHZACC(1) - PHTHRS(1) ! Excess thermal time
PHZACC(2) = ExcsTTD ! Initializes Phase 2
ExcsTTD = 0 ! Voids excess thermal time
END IF
ELSE IF (CropTypeCode >= 2) THEN
PHZACC(1) = PHZACC(1) + DeltaTTD
IF (PHZACC(1) > PHTHRS(1)) THEN
! Stage 2 is skipped
PhenoStage = 3 ! Stage 3 is set right away
!-----------------------------------------------------------------------
! Emergence, next stage, occurs on day DAS
!-----------------------------------------------------------------------
NVEG0 = DAS
YREMRG = YRDOY
STGDOY(2) = YRDOY
STGDOY(3) = YRDOY
DayOfStage(2) = YRDOY ! Phases 2 and 3 occur the same day
DayOfStage(3) = YRDOY ! Phases 2 and 3 occur the same day
!-----------------------------------------------------------------------
DaysAfterPlantOfStage(2) = DAP ! Idem
DaysAfterPlantOfStage(3) = DAP ! Idem
ExcsTTD = PHZACC(1) - PHTHRS(1) ! Excess thermal time
! Initializes thermal time accumulation for Phase 3
PHZACC(3) = ExcsTTD
ExcsTTD = 0 ! Voids excess thermal time
PHZACC(2) = 0 ! Second phase does not take place
! Primary stalk number for ratoon depends on mature stalks harvested
! during the previous season.
! Starts development of stalks; NewStalk becomes 1 now
NewStalk = StkHrNO * RTNFAC
DO Stalk = 1, NewStalk
StalkState(Stalk,1) = 'LIVE'
StalkState(Stalk,2) = 'PRIM'
END DO
NewStalk = StkHrNO * RTNFAC !(temp for debugging)
END IF
END IF
!=============
CASE (2)
!=============
IF (CropTypeCode == 1) THEN
PHZACC(2) = PHZACC(2) + DeltaTTD
DepthToEmer = DepthToEmer + DeltaDepthOfEmer
IF (DepthToEmer > SDEPTH) THEN
PhenoStage = 3
!-----------------------------------------------------------------------
! Emergence, next stage, occurs on day DAS
!-----------------------------------------------------------------------
NVEG0 = DAS
YREMRG = YRDOY
STGDOY(3) = YRDOY
DayOfStage(3) = YRDOY
!-----------------------------------------------------------------------
DaysAfterPlantOfStage(3) = DAP
ExcsTTD = PHZACC(2) - PHTHRS(2) ! Excess thermal time
! Initializes thermal time accumulation for Phase 3
PHZACC(3) = ExcsTTD
ExcsTTD = 0 ! Voids excess thermal time
! Starts development of stalks; NewStalk must become 1 now
NewStalk = NewStalk + 1
StalkState(NewStalk,1) = 'LIVE'
StalkState(NewStalk,2) = 'PRIM'
!Only 1st stalk is counted here
END IF
!! ELSE IF (CropTypeCode >= 2) THEN ! Ratoon crop
!! PHZACC(2) = 0 ! Second phase does not take place
!! PhenoStage = 3
!-----------------------------------------------------------------------
! Emergence, next stage, occurs on day DAS
!-----------------------------------------------------------------------
!! NVEG0 = DAS
!! YREMRG = YRDOY
!! STGDOY(3) = YRDOY
!! DayOfStage(3) = YRDOY
!-----------------------------------------------------------------------
! Primary stalk number for ratoon depends on mature stalks harvested
! during the previous season.
!! DaysAfterPlantOfStage(3) = DAP
!************************************************************************
! Check this section: is it ever used?
! Starts development of stalks; NewStalk must become 1 now
!! NewStalk = NewStalk + 1
!! StalkState(NewStalk,1) = 'LIVE'
END IF
!************************************************************************
!=============
CASE (3)
!=============
PHZACC(3) = PHZACC(3) + DeltaTTD
!************* Stalk-based stalk development code 3(FSR) ****************
DO Stalk = 1, Smax
IF (StalkState(Stalk,1) .EQ. 'LIVE') THEN
NewTiller = NewTiller
& + DeltaTillerNum(Stalk)
& + TillerExcess
TillerExcess = 0.
ELSE
END IF ! (StalkState(Stalk,1) .EQ. 'LIVE')
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
IF (INT(NewTiller) >= 1.0) THEN
NewStalk = NewStalk + INT(NewTiller)
StalkState(NewStalk,1) = 'LIVE'
StalkState(NewStalk,2) = 'TILR'
TillerExcess = NewTiller - INT(NewTiller)
TillerCount = Tillercount + 1
NewTiller = 0
END IF ! ((INT(NewTiller) >= 1.0)
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Temporary limit
IF (NewStalk > Smax) THEN
StalkState(NewStalk,1) = ' '
NewStalk = Smax
END IF
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Remove stalks not receiving sufficient CH2O for maint respiration
IF (Kill(Stalk) .EQ. 1) THEN
StalkState(Stalk,1) = 'DEAD'
Kill(Stalk) = -1
EXIT ! Exits DO Loop to avoid senescing more than one stalk
ELSE ! per day
END IF
END DO
!************* Stalk-based stalk development code 3(FSR) ****************
! Section for leaf development
! Integrates thermal time for each day and each stalk that has been set
DO Stalk = 1, Smax
IF (StalkState(Stalk,1) .EQ. 'LIVE') THEN
SumTTStalk(DAS, Stalk) =
& SumTTStalk(DAS - 1, Stalk) + DeltaTTD
IF (LeafNum(DAS-1, Stalk) .EQ. 0) THEN !to start stalk with
LeafNum(DAS-1, Stalk) = 1 ! a leaf
END IF
ELSE
SumTTStalk(DAS, Stalk) = 0
END IF
! Integrates number of leaves in each day and each stalk
LeafNum(DAS, Stalk) =
& LeafNum(DAS - 1, Stalk) + DeltaLeafNum(DAS, Stalk)
END DO
IF (PHZACC(3) > PHTHRS(3)) THEN
PhenoStage = 4
STGDOY(4) = YRDOY
DayOfStage(4) = YRDOY
DaysAfterPlantOfStage(4) = DAP
ExcsTTD = PHZACC(3) - PHTHRS(3) ! Excess thermal time
! Initializes thermal time accumulation for Phase 3
PHZACC(4) = ExcsTTD
ExcsTTD = 0 ! Voids excess thermal time
END IF
!=============
CASE (4)
!=============
PHZACC(4) = PHZACC(4) + DeltaTTD
!************* Stalk-based stalk development code 4(FSR) ****************
DO Stalk = 1, Smax
IF (StalkState(Stalk,1) .EQ. 'LIVE') THEN
NewTiller = NewTiller
& + DeltaTillerNum(Stalk)
& + TillerExcess
TillerExcess = 0.
ELSE
END IF ! (StalkState(Stalk,1) .EQ. 'LIVE')
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
IF (INT(NewTiller) >= 1.0) THEN
NewStalk = NewStalk + INT(NewTiller)
StalkState(NewStalk,1) = 'LIVE'
StalkState(NewStalk,2) = 'TILR'
TillerExcess = NewTiller - INT(NewTiller)
TillerCount = Tillercount + 1
NewTiller = 0
END IF ! ((INT(NewTiller) >= 1.0)
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Temporary limit
IF (NewStalk > Smax) THEN
StalkState(NewStalk,1) = ' '
NewStalk = Smax
END IF
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Remove stalks not receiving sufficient CH2O for maint respiration
IF (Kill(Stalk) .EQ. 1) THEN
StalkState(Stalk,1) = 'DEAD'
Kill(Stalk) = -1
EXIT ! Exits DO Loop to avoid senescing more than one stalk
ELSE ! per day
END IF
END DO
! - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Check for live stalks. If none are found, season (run) ends.
livecount = 0
DO Stalk = 1,Smax
IF (StalkState(Stalk,1) .EQ. 'LIVE') THEN
livecount = livecount +1
END IF !(StalkState(Stalk,1) .EQ. 'LIVE') etc
END DO ! Stalk
IF (livecount .LT. 1) THEN
!!! MDATE = YRDOY ! punch line goes here
livecount = livecount ! temp for debugging
END IF !
!************* Stalk-based stalk development code 4(FSR) ****************