-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathreapar.f90
1767 lines (1554 loc) · 72.6 KB
/
reapar.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
SUBROUTINE REAPAR(CBINFLAG,AEFFA,BETAMI,BETAMX,CALPHA,CFLEPS,CFLPAR, &
CFLSHP,CMDSOL,CMDFFT,CMDFRM,CSHAPE,CMDTRQ,CXE01_LF, &
CXE02_LF,DEGRAD,ETASCA,EXTNDXYZ,GAMMA,IDVERR,IDVOUT, &
INIT,IOPAR,IORTH,JPBC,IWRKSC,IWRPOL,MXITER,MXNX,MXNY, &
MXNZ,MXBETA,MXCOMP,MXPHI,MXRAD,MXSCA,MXTHET,MXWAV, &
NBETA,NCOMP,NPHI,NRAD,NRFLD,NRFLDB,NSCA,NTHETA,NWAV, &
IWAV0,IRAD0,IORI0,NSMELTS,PHIMIN,PHIMAX,PHIN,THETAN, &
THETMI,THETMX,TOL,WAVEA,SHPAR,SMIND1,SMIND2,DX,NAMBIENT,CENTER,ELENERGY)
!Argument 'CENTER' and 'ELENERGY' added by SMC 8.5.13 following NWB 3/13/12
USE DDPRECISION,ONLY: WP
USE DDCOMMON_0,ONLY: IDIPINT
IMPLICIT NONE
!----------------------------- reapar_v3 -----------------------------------
! Arguments:
CHARACTER :: CALPHA*(*),CBINFLAG*(*),CFLPAR*(*),CFLSHP*80,CMDSOL*(*), &
CMDFFT*(*),CMDFRM*(*),CMDTRQ*(*),CSHAPE*(*)
INTEGER :: IDVERR,IDVOUT,INIT,IOPAR,IORI0,IORTH,JPBC,IRAD0,IWAV0, &
IWRKSC,IWRPOL,MXBETA,MXCOMP,MXITER,MXNX,MXNY,MXNZ,MXPHI, &
MXRAD,MXSCA,MXTHET,MXWAV,NBETA,NCOMP,NRFLD,NRFLDB,NPHI, &
NRAD,NSCA,NSMELTS,NWAV,NTHETA
INTEGER :: &
SMIND1(9), &
SMIND2(9)
CHARACTER(60) :: &
CFLEPS(MXCOMP)
REAL(WP) :: BETAMI,BETAMX,DEGRAD,ETASCA,GAMMA,NAMBIENT,PHIMIN,PHIMAX, &
THETMI,THETMX,TOL,ELENERGY !ELENERGY added by SMC 8.5.13 following NWB 7/12/12
REAL(WP) :: &
AEFFA(1:MXRAD), &
DX(1:3), &
EXTNDXYZ(1:6), &
PHIN(MXSCA), &
SHPAR(1:12), &
THETAN(MXSCA), &
WAVEA(MXWAV), &
CENTER(3) !CENTER added by SMC 8.5.13 following NWB 3/13/12
COMPLEX(WP) :: &
CXE01_LF(1:3), &
CXE02_LF(1:3)
! Local variables:
CHARACTER :: CDIVID*3,CLINE*70,CMSGNM*70,CWHERE*70
INTEGER :: J,JPLANE,NF,NPLANES,NSCA0
REAL(WP) :: AEFEND,AEFINI,DELTA,DTHETA,E1,PHI1,THETA1,THETA2, &
WAVEND,WAVINI
! External routines:
EXTERNAL DIVIDE,ERRMSG,WRIMSG
!***********************************************************************
! Subroutine REAPAR handles the reading of input parameters from the
! "ddscat.par" file, as well as elementary processing with those input
! parameters to generate arrays.
! Given:
! CFLPAR = name of file (normally 'ddscat.par')
! NRFLD = 0 on first call
! = 1 if previous call set NRFLD=1 to prepare to do
! nearfield calculation
! Returns:
! CMDTRQ = 'NOTORQ' or 'DOTORQ'
! CMDSOL = 'PBCGST' or 'PBCGS2' or 'PETRKP' or 'GPBICG' or 'QMRCCG' or
! 'SBIGCG'
! CMDFFT = 'GPFAFT' or 'FFTW21' or 'FFTMKL'
! CALPHA = 'LATTDR' or 'GKDLDR' or 'FLTRCD'
! CBINFLAG = 'ALLBIN' or 'ORIBIN' or 'NOTBIN'
! CSHAPE = 'RCTGLPRSM' or 'ELLIPSOID' or ... other valid shape
! JPBC = 0 :isolated finite targets
! = 1 :target periodic in y direction, finite in z
! = 2 :target finite in y direction, periodic in z
! = 3 :target periodic in both y and z directions
! NCOMP = number of compositions
! CFLEPS(1-NCOMP) = filenames for dielectric functions
! NRFLD = 0 : skip near-field calculations
! = 1 : to prepare to do additional near-field calculation
! = 2 : to perform near-field calculation itself
! IDIPINT = 0 for point dipole interaction method
! = 1 for filtered coupled dipole (FCD) interaction method
! NRFLDB = 0 : to skip near-field calculation of B
! = 1 : to perform near-field calculation of B
! EXTNDXYZ(1-6) = fractional expansion of calculational volume
! in -x,+x,-y,+y,-z,+z directions
! (used only if NRFLD=1)
! INIT = 0 (no longer used)
! TOL = error tolerance
! GAMMA = parameter controlling integration limit for PBC
! ETASCA = parameter controlling accuracy of calculation
! of <cos(theta)>
! 1 is OK, 0.5 gives high accuracy
! NWAV = number of wavelengths
! WAVEA(1-NWAV) = wavelengths (in vacuo, physical units)
! NAMBIENT = (real) refractive index of ambient medium
! NRAD = number of radii
! AEFFA(1-NRAD) = target effective radii (physical units)
! CXE01_LF(1-3)= (complex) polarization state 1
! CXE02_LF(1-3)= complex polarization state 2 (orthogonal to 1)
! IORTH = 1 or 2 (number of incident polarization states to use
! IWRKSC = 0 or 1 (not write/write ".sca" files)
! IWRPOL = 0 or 1 (not write/write ".pol" files)
! NBETA = number of beta values for target rotation
! BETAMI = minimum beta value (rad) [input from file in deg]
! BETAMX = maximum beta value (rad) [input from file in deg]
! NTHETA = number of theta values for target rotation
! THETAMI = minimum theta value (rad) [input from file in deg]
! THETAMX = maximum theta value (rad) [input from file in deg]
! NPHI = number of PHI values for target rotation
! PHIMIN = minimum PHI value (rad) [input from file in deg]
! PHIMAX = maximum PHI value (rad) [input from file in deg]
! NSMELTS = number of elements of S matrix to calculate
! (less than or equal to 9)
! SMINDI(J) = indices of S matrix elements to be calculated
! (e.g., 11 , 21 , 31 , 41 , 12 , 13)
! CMDFRM = 'LFRAME' or 'TFRAME' (Lab frame or Target frame)
! NSCA = number of scattering directions
! If JPBC=0:
! THETAN(1-NSCA) = theta for scattering directions (rad)
! PHIN(1-NSCA) = phi for scattering directions (rad)
! If JPBC=1 or 2:
! THETAN(1-NSCA) = order_y or order_z for scattering directions
! PHIN(1-NSCA) = psi (radians) for scattering directions,
! where psi = rotation around axis of
! target periodicity (y_TF if JPBC=1, z_TF if JPB
! If JPBC=3: periodic in both y and z directions
! THETAN(1-NSCA) = order_y for diffraction
! PHIN(1-NSCA) = order_z for diffraction
! Original version created by P.J.Flatau, Colorado State Univ.
! Modified by B.T.Draine, Princeton Univ. Obs.
! History:
! 90.11.09 (BTD): Modified so as to remove sensitivity to errors
! in sign of entered value of DTHETA.
! 90.11.30 (BTD): Added test to compare CMDFFT and MXMEM
! 90.12.14 (BTD): Changed option MKRECT to RCTGLPRSM
! 91.01.03 (BTD): Added MXBETA,MXPHI,MXTHET to argument list
! and added checking of NBETA,NPHI,NTHETA
! 91.05.08 (BTD): Added CALPHA to allow choice of prescription
! for polarizabilities
! 91.05.08 (BTD): Some rewriting done.
! 91.05.14 (BTD): Added separate options LDRXYZ and LDRAVG
! 91.05.15 (BTD): Added new option LDR000
! 91.05.23 (BTD): Added variable IWRKSC to argument list and
! added code to read IWRKSC
! 91.05.30 (BTD): Changed LDRXYZ->LDR100
! Added option LDR111
! 91.11.13 (BTD): Now three options: DRAI88, GOBR88, LATTDR
! (subroutine ALPHA now allows for directional
! dependence for LATTDR option)
! 92.05.14 (BTD): Allow for experimental option LDRISO
! 92.09.21 (BTD): Allow for option ANIELLIPS (ellipsoid of anisotropic
! material); change option GRPHIT -> UNIAXICYL (cylinder
! of uniaxial material)
! 93.01.07 (BTD): Allow for options ELLIPSO_2 and ANI_ELL_2 (two touching
! ellipsoids, each either isotropic or anisotropic)
! 93.03.11 (BTD): Deleted all code associated with variable CMACHN
! which previously was used to identify machine/OS
! This information was not used or needed.
! 93.12.15 (BTD): Added code to check for compatibility of NCOMP
! when CSHAPE=ELLIPSO_2,ELLIPSO_3,ANI_ELL_2,or ANI_ELL_3
! 94.01.27 (BTD): Replaced SHPAR1,SHPAR2,SHPAR3 by SHPAR(1-6)
! Added CONELLIPS (two concentric ellipsoids) as target
! option.
! 94.06.20 (PJF): Added code to handle NEWTMP
! 94.12.02 (BTD): Added options
! LB1991 (Dungey-Bohren 1991 polarizability)
! LR1994 (Lumme-Rahola 1994 polarizability)
! 95.05.26 (BTD): Corrected construction of orthogonal polarization
! state for case where E01 is complex (e.g.,
! circular polarization)
! 95.06.14 (BTD): Added code to check that user is not requesting
! f_ml when illuminating target with other then
! linearly polarized light
! 95.06.19 (BTD): Added argument CMDSOL to allow selection of
! different CCG methods
! Added argument CMDTRQ to allow torque calculations
! to be skipped.
! 96.01.05 (BTD): Added shape options MLTBLOCKS and DW1996TAR
! 96.01.25 (BTD): Added shape option SPHROID_2
! 96.01.25 (BTD): Disable code added 95.06.14 so that user can
! calculate f_ml even when illuminating target with
! elliptically polarized light.
! 96.10.18 (BTD): Changed option NEWTMP to GPFAFT (Generalized Prime
! Factor Algorithm for Fourier Transform).
! 96.11.06 (BTD): Added provision to
! (1) produce fatal error if CXE01(1).NE.(0.,0.)
! (2) normalize CXE01 if necessary
! 96.11.14 (PJF): cbinflag for binary file
! 96.11.15 (PJF): netCDF netflag
! 96.11.21 (BTD): modified handling of CBINFLAG and CNETFLAG to
! be consistent with other option flags
! 97.11.01 (BTD): add DX to argument list
! read DX from ddscat.par
! renormalize DX (input DX not required to be
! normalized with DX(1)*DX(2)*DX(3)=1.)
! 98.12.06 (BTD): modified to read NSMELTS,SMELTS from ddscat.par
! and pass NSMELTS,SMIND1,SMIND2 as arguments to calling
! program DDSCAT
! 98.12.21 (BTD): changed dimension of CFLPAR from CHARACTER*40 to
! CHARACTER*60 to allow longer file names
! (also changed in DDSCAT.f and dielec.f)
! 98.12.29 (BTD): modified to support new option LYRD_SLAB
! 99.01.26 (BTD): modified to allow first IWAV,IRAD,IORI to be
! specified in ddscat.par (useful for restarting
! jobs which did not complete)
! 99.06.29 (BTD): corrected input of SHPAR for option SPHROID_2
! previously had input only SHPAR(1-5),
! now input SHPAR(1-6)
! 00.06.12 (BTD): modified to support SPHERES_N option
! 00.06.23 (BTD): modified to support FFTWFJ option
! remove MXMEM from argument list
! since no longer supporting TMPRTN FFT option which
! required additional memory
! 00.07.05 (BTD): modified to print more information about CMDFFT
! option
! 02.02.12 (BTD): modified to support TRNGLPRSM option
! 03.01.29 (BTD): modified to force DX(1-3)=1.
! 03.05.21 (BTD): modified to support SPHARM option
! 03.05.22 (BTD): bug fix re SPHARM option
! 03.07.13 (BTD): changed FFTWFJ to FFTW21 in anticipation of future
! option FFTW3x (FFTW 2.1.x and FFTW 3.0.1 have
! different interfaces, so they must be distinguished).
! 03.10.23 (BTD): removed ICTHM and IPHIM from argument list,
! removed code to read these values from ddscat.par
! 03.10.24 (BTD): Add parameter ETASCA to allow user to control
! angular resolution used for scattering calculations.
! Add checks against invalid input
! 04.02.26 (BTD): Modify to accept polarizability option 'GKDLDR'
! for Gutkowicz-Krusin & Draine (2004) lattice
! dispersion relation.
! 04.04.01 (BTD): Modify to support option NPSPBC (multi-sphere cluster
! with periodic boundary conditions)
! 04.04.05 (BTD): Remove SPHARM option: GAUSS_SPH is sufficient
! 04.04.29 (BTD): Added ANIRCTNGL option for homogeneous anisotropic
! rectangular target
! 04.05.23 (BTD): Added CWHERE string to provide error information
! in case of malformed ddscat.par file.
! 04.09.14 (BTD): Permit new target option ANIFRMFIL, to read target data
! from file, with target allowed to be anistropic
! with arbitrary orientations of local "Dielectric Frame
! (in which local dielectric tensor is diagonalized)
! relative to Target Frame.
! 04.10.12 (BTD): Modified to support new target option SPH_ANI_N, for
! a target composed of N possibly anisotropic spheres.
! 04.10.14 (BTD): Modified to support option of reading wavelengths
! from file 'wave.tab' if CDIVID='TAB' in ddscat.par
! 05.03.19 (BTD): Added 2 lines to set CLFEPS(1) in case using options
! H2OLIQ or H2OICE (previously were left undefined,
! which cause problems with output statement)
! 05.08.04 (BTD): Added requirement to read string CMDFRM where
! CMDFRM = 'LFRAME' to specify scattering directions in
! Lab Frame (rel. to xlab,ylab,zlab,
! where xlab = direction of incident beam
! CMDFRM = 'TFRAME' to specify scattering directions in
! Target Frame (rel. to a1,a2,a3)
! 05.10.18 (BTD): Added check to ensure that NBETA*NTHETA*NPHI < 1001
! if IWRKSC=1 (since current version of NAMER is limited
! to maximum of 1000 directions)
! 06.04.10 (BTD): Added IWRPOL to argument list, and modified REAPAR
! to read IWRPOL from ddscat.par in order to allow
! user to elect to write or to not write out polarizatio
! array for subsequent processing.
! 06.09.13 (BTD): Added handling of option CYLNDRCAP
! 06.09.15 (BTD): Modified to support new version of target routine
! TARHEX, that now requires SHPAR(3) to select prism
! orientation.
! 06.09.15 (BTD): Modified to support new target option HEXGONPBC
! 06.10.04 (BTD): Extension to support new scattering directions
! required for far-field scattering when PBC is used.
! 06.10.11 (BTD): Corrected errors in handling of JPBC=2,3,
! 06.11.29 (BTD): Added code to catch error in input CMDFRM
! 06.12.08 (BTD): Added support for target option DSKRCTPBC
! 07.01.20 (BTD): Added support for target option DSKRCTNGL
! 07.02.23 (BTD): Added support for target option LYRSLBPBC
! 07.02.25 (BTD): Corrected bug (typo)
! 07.08.05 (BTD): Version 7.0.3
! Add code to read MXNX,MXNY,MXNZ from ddscat.par
! 07.08.12 (BTD): Allow CMDSOL='PBCGS2'
! 07.10.26 (BTD): Add support for target option DSKBLYPBC
! reorder input shape paramers SHPAR for DSKRCTPBC
! changed SHPAR(6) -> SHPAR(10)
! 07.10.27 (BTD): Further improve of reapar
! Eliminate CDIEL -- option 'TABLES' is now the only
! one used, so no need to input and carry around CDIEL
! 08.02.01 (BTD): Changed SHPAR(10) -> SHPAR(12) to allow up to 12
! shape parameters to be input
! Added support for target option BLSLINPBC
! Added support for target option RCTGLBLK3
! Added support for target option TRILYRPBC
! 08.03.11 (BTD): corrected typo: BLSLINPBC -> BISLINPBC
! 08.03.11 (BTD): ver7.0.5
! Added ALPHA to argument list.
! Modified to read parameter ALPHAA from ddscat.par
! 08.03.12 (BTD): Modified to read NPLANES from ddscat.par
! so that we no longer rely on END condition when reading
! scattering directions or scattering planes from ddscat.par
! 08.03.12 (BTD): Modified to use new form of ddscat.par
! eliminate reading of INIT (since we always set INIT=0)
! add new descriptive lines separating ddscat.par into sections
! 08.04.19 (BTD): Changed notation: ALPHA -> GAMMA
! reordered argument list
! 08.06.07 (BTD): Added support for FFT option FFTMKL
! 08.07.22 (BTD): Added sanity check to catch users who call for SPHERES_N but
! specify multiple dielectric functions
! 08.08.29 (BTD): Removed CNETFLAG from argument list
! Removed code to read CNETFLAG from ddscat.par
! 08.09.12 (BTD): Corrected typo DSKRCTGNL -> DSKRCTNGL
! 08.09.16 (BTD): Corrected bug reported by Wang Lin, Nankai Univ.:
! typo SPHERN_PBC -> SPHRN_PBC caused SPHERN_PBC to be
! assigned JPBC=0 regardless of values of SHPAR(2) or SHPAR(3).
! 08.09.23 (BTD): Corrected typo DSKBLYPRC -> DSKLBYLRPBC line 613
! [typo reported by Georges Levi, Univ-Paris-Diderot]
! 09.08.27 (BTD): Added code to test for insane value of TOL
! (which could result if other than NCOMP diel. fns.
! are listed in ddscat.par)
! 09.09.10 (BTD): ver7.1.0
! Added support for new target options
! ANIFILPBC
! FRMFILPBC
! Added ability to recognize lowercase target options
! 09.12.11 (BTD): Corrected coding error which prevented use of
! target options DW1996TAR and MLTBLOCKS
! 10.01.02 (BTD): Add code to write out SHPAR info to running output
! 10.02.05 (BTD): Corrected error that caused DW1996TAR option to fail
! 10.05.08 (BTD): ver7.2.0
! * Added support for CCG option GPBICG
! * Added support for CCG option QMRCCG
! 10.05.09 (BTD): * set MXITER in ddscat.par
! 11.07.29 (BTD): * Add support for target option NEARFIELD
! * corrected error affecting reading input for SPHROID_2
! 11.08.12 (BTD): ver7.2.1
! * eliminate target option NEARFIELD
! * Add NRFLD and EXTNDXYZ to argument list
! * Add code to read NRFLD and EXTNDXYZ
! 11.08.16 (BTD): * suppress diagnostic write statements
! 11.08.30 (BTD): ver7.2.2
! * add support for refractive index of ambient medium
! * add variable NAMBIENT to argument list
! 11.08.31 (BTD) * improve diagnostics
! 11.10.18 (BTD) * add support for target optino EL_IN_RCT
! 11.11.07 (BTD) ver7.1.1 and ver7.2.2
! * Corrected error reported by V. Choliy 11.11.07
! had failed convert PHI1 from degrees to radians
! for JPBC=1 or JPBC=2
! 11.11.18 (BTD) * change notation
! CXE01 -> CXE01_LF
! CXE02 -> CXE02_LF
! * added support for reading aeff from table
! 12.01.23 (BTD) * suppress line to read IWRPOL from ddscat.par
! * permanently set IWRPOL=0
! 12.07.10 (BTD) reapar_v2 and ver7.3
! * added NRFLDB to argument list
! * added code to set NRFLDB=0 if NF=0 or 1 in ddscat.par
! set NRFLDB=1 if NF=2
! 12.08.02 (IYW) * added DIPINT to argument list
! 12.08.11 (BTD) * removed DIPINT from argument list
! * added IDIPINT from module DDCOMMON_0
! 12.12.28 (BTD) reapar_v3
! * eliminate reading IDIPINT from ddscat.par
! * added option to read DDA method = FLTRCD
! and set IDIPINT=1
! 13.01.04 (BTD) * corrected typo aeffa(nwav) -> aeffa(nrad)
! 13.01.10 (BTD) * added CG option SBICGM
! end history
! Copyright (C) 1993,1994,1995,1996,1997,1998,1999,2000,2002,2003,2004
! 2005,2006,2007,2008,2009,2010,2011,2012,2013
! B.T. Draine and P.J. Flatau
! This code is covered by the GNU General Public License.
!***********************************************************************
! 110801 (BTD) variable INIT no longer used, but definite it anyway
INIT=0
OPEN(UNIT=IOPAR,FILE=CFLPAR,STATUS='OLD')
CWHERE='error reading line 1 of ddscat.par'
READ(IOPAR,FMT=*,ERR=99)CLINE
CALL WRIMSG('REAPAR',CLINE)
CWHERE='error reading line 2 of ddscat.par'
READ(IOPAR,FMT=*,ERR=99)CLINE
CALL WRIMSG(' ',CLINE)
!=======================================================================
! Specify whether torques are to be calculated:
! CMDTRQ*6 = 'DOTORQ' - calculate torque on grain
! = 'NOTORQ' - do not calculate torque on grain
CWHERE='error reading CMDTRQ from ddscat.par'
READ(IOPAR,FMT=*,ERR=99)CMDTRQ
IF(CMDTRQ=='dotorq')CMDTRQ='DOTORQ'
IF(CMDTRQ=='notorq')CMDTRQ='NOTORQ'
IF(CMDTRQ=='DOTORQ')THEN
WRITE(CMSGNM,FMT='(A,A)')CMDTRQ,' - compute torques '
CALL WRIMSG('REAPAR',CMSGNM)
ELSEIF(CMDTRQ=='NOTORQ')THEN
WRITE(CMSGNM,FMT='(A,A)')CMDTRQ,' - do not compute torques '
CALL WRIMSG('REAPAR',CMSGNM)
ELSE
!*** diagnostic
! write(0,*)'reapar ckpt 0.1 : CMDTRQ=',CMDTRQ
!***
CALL ERRMSG('FATAL','REAPAR',' wrong definition of CMDTRQ')
ENDIF
!=======================================================================
! Define method used for iterative solution of complex linear equations
! CMDSOL*6 = 'PETRKP' - Petravic & Kuo-Petravic method
! = 'PBCGST' - PIM BiConjugate Gradient with Stabilization
! = 'PBCGS2' - M.A Botcheve implementation of BiCGstab
! enhanced to improve convergence properties
! with finite precision arithmetic
! = 'GPBICG' - Tang et al conjugate gradient solver implemented
! by P.C. Chaumet & A. Rahmani
! = 'QMRCCG' - PIM interface to QMR solver implemented by
! P.C. Chaumet and A. Rahmani
! = 'SBICGM' - Simplified Bi Conjugate Gradient Method
! implemented by Sarkar, Yang, & Arvas 1988
CWHERE='error reading CMDSOL from ddscat.par'
READ(IOPAR,FMT=*,ERR=99) CMDSOL
IF(CMDSOL=='petrkp')CMDSOL='PETRKP'
IF(CMDSOL=='pbcgst')CMDSOL='PBCGST'
IF(CMDSOL=='pbcgs2')CMDSOL='PBCGS2'
IF(CMDSOL=='gpbicg')CMDSOL='GPBICG'
IF(CMDSOL=='qmrccg')CMDSOL='QMRCCG'
IF(CMDSOL=='sbicgm')CMDSOL='SBICGM'
IF(CMDSOL=='PETRKP'.OR. &
CMDSOL=='PBCGST'.OR. &
CMDSOL=='PBCGS2'.OR. &
CMDSOL=='GPBICG'.OR. &
CMDSOL=='QMRCCG'.OR. &
CMDSOL=='SBICGM')THEN
WRITE(CMSGNM,FMT='(A,A)')CMDSOL,' - CCG Method '
CALL WRIMSG('REAPAR',CMSGNM)
ELSE
CALL ERRMSG('FATAL','REAPAR',' wrong definition of CMDSOL')
ENDIF
!=======================================================================
! Define FFT method:
! CMDFFT*6 = 'GPFAFT' - GPFA code of Temperton
! = 'FFTW21' - FFTW 2.1.x code of Frigo & Johnson
! = 'FFTMKL' - Use DFTI from Intel MKL
CWHERE='error reading CMDFFT from ddscat.par'
READ(IOPAR,FMT=*,ERR=99) CMDFFT
IF(CMDFFT=='FFTW21')THEN
WRITE(CMSGNM,FMT='(A,A,A)')CMDFFT,' - using FFTW 2.1.x ', &
'package from Frigo & Johnson'
CALL WRIMSG('REAPAR',CMSGNM)
ELSEIF(CMDFFT=='GPFAFT')THEN
WRITE(CMSGNM,FMT='(A,A,A)')CMDFFT,' - using GPFA package ', &
'from Clive Temperton'
CALL WRIMSG('REAPAR',CMSGNM)
ELSEIF(CMDFFT=='FFTMKL')THEN
WRITE(CMSGNM,FMT='(A,A,A)')CMDFFT,' - using DFTI from Intel ', &
'Math Kernel Library (MKL)'
CALL WRIMSG('REAPAR',CMSGNM)
ELSE
WRITE(CMSGNM,FMT='(A,A)')'DDSCAT 7.0 only supports FFT options', &
' FFTW21, GPFAFT, and FFTMKL'
CALL WRIMSG('REAPAR',CMSGNM)
CALL ERRMSG('FATAL','REAPAR',' wrong definition of CMDFFT')
ENDIF
!=======================================================================
! Define prescription for computing polarizabilities:
! CALPHA*6 = 'LATTDR' - Lattice Dispersion Relation of
! Draine & Goodman (1993)
! 'GKDLDR' - Lattice Dispersion Relation of
! Gutkowicz-Krusin & Draine (2004)
! 'FLTRCD' - Filtered Discrete Dipole treatment of
! Piller & Martin (1998) as corrected by
! Gay-Balmaz & Martin (2002) and discussed by
! Yurkin, Min, & Hoekstra (2010)
CWHERE='error reading CALPHA from ddscat.par'
READ(IOPAR,FMT=*,ERR=99)CALPHA
IF(CALPHA=='LATTDR')THEN
CALL WRIMSG('REAPAR','LATTDR - Draine & Goodman (1993) LDR for alpha')
IDIPINT=0
ELSEIF(CALPHA=='GKDLDR')THEN
CALL WRIMSG('REAPAR', &
'GKDLDR - Gutkowicz-Krusin & Draine (2004) LDR for alpha')
IDIPINT=0
ELSEIF(CALPHA=='FLTRCD')THEN
CALL WRIMSG('REAPAR', &
'FLTRCD - Filtered Coupled Dipole (Piller & Martin 1998)')
IDIPINT=1
ELSE
CALL ERRMSG('FATAL','REAPAR',' wrong definition of CALPHA')
ENDIF
!=======================================================================
! binary file flag
CWHERE='error reading CBINFLAG from ddscat.par'
READ(IOPAR,FMT=*,ERR=99)CBINFLAG
IF(CBINFLAG=='allbin')CBINFLAG='ALLBIN'
IF(CBINFLAG=='oribin')CBINFLAG='ORIBIN'
IF(CBINFLAG=='notbin')CBINFLAG='NOTBIN'
IF(CBINFLAG=='ALLBIN'.OR.CBINFLAG=='ORIBIN'.OR.CBINFLAG=='NOTBIN')THEN
WRITE(CMSGNM,FMT='(A,A)')CBINFLAG,' - Unformatted binary dump option'
CALL WRIMSG('REAPAR',CMSGNM)
ELSE
CALL ERRMSG('FATAL','REAPAR',' Wrong definition of CBINFLAG')
ENDIF
!=======================================================================
! Read upper bound on target extent
! skip line:
READ(IOPAR,FMT=*,ERR=99)CLINE
CALL WRIMSG('REAPAR',CLINE)
CWHERE='error reading target size bounds MXNX,MXNY,MXNZ from ddscat.par'
READ(IOPAR,FMT=*,ERR=99)MXNX,MXNY,MXNZ
WRITE(CMSGNM,FMT='(A,3I6,A)')'allow MXNX,MXNY,MXNZ=',MXNX,MXNY,MXNZ, &
' for target generation'
CALL WRIMSG('REAPAR',CMSGNM)
!=======================================================================
! e-BEAM CENTROID READER
! Added by NWB 3/13/12
READ(IOPAR,FMT=*,ERR=99)CENTER
!=======================================================================
! ELECTRON RELATIVISTIC ENERGY READER
! Added by SMC8.5.13 following NWB 7/12/12
READ(IOPAR,FMT=*,ERR=99)ELENERGY
!=======================================================================
! Define shape:
! skip line:
READ(IOPAR,FMT=*,ERR=99)CLINE
CALL WRIMSG('REAPAR',CLINE)
! CSHAPE*9 = 'FROM_FILE' shape and composition data will later
! be read from file CFLSHP
! for dielectric tensors that are diagonal
! in Target Frame
! = 'ANIFRMFIL' read shape and composition data from file,
! for general anisotropic dielectric tensors
! with orientation angles THETADF,PHIDF,BETAD
! relative to Target Frame.
! = 'FRMFILPBC' shape and composition data for TUC will later
! be read from file CFLSHP
! for dielectric tensors that are diagonal
! in Target Frame
! PYD and PZD are input via ddscat.par
! = 'ANIFILPBC' shape and composition data for TUC will later
! be read from file CFLSHP
! for general anisotropic dielectricc tensors
! with orientation angles THETADF,PHIDF,BETAD
! relative to Target Frame.
! PYD and PZD are input via ddscat.par
! = 'ANIELLIPS' ellipsoid of anisotropic material
! = 'ANIRCTNGL' homogeneous anisotropic rectangular target
! = 'ANI_ELL_2' two touching anisotropic ellipsoids of
! materials 1-6
! = 'ANI_ELL_3' three touching anisotropic ellipsoids of
! materials 1-9
! = 'BISLINPBC' bilayer slab with periodic grid of lines
! parallel to z on top, with y-period/d=PYD
! [or, if PYD=0, a single line]
! = 'CONELLIPS' two concentric ellipsoids of materials 1,2
! = 'CYLINDER1' homogeneous finite cylinder
! = 'CYLNDRCAP' homogeneous cylinder with hemispherical endcaps
! = 'CYLNDRPBC' 1-d or 2-d array of finite cylinders
! = 'DSKBLYPBC' 1-d or 2-d array of disk on bilayer rect. slab
! = 'DSKRCTNGL' single disk on rectangular slab
! = 'DSKRCTPBC' 1-d or 2-d array of disk on rectangular slab
! = 'DW1996TAR' 13-cube target used by Draine & Weingartner 1996
! = 'EL_IN_RCT' ellipsoid embedded in rectangular block
! = 'ELLIPSOID' ellipsoid (homogeneous and isotropic)
! = 'ELLIPSPBC' 1-d or 2-d array of ellipsoids
! = 'ELLIPSO_2' two touching isotropic ellipsoids of
! materials 1 and 2
! = 'ELLIPSO_3' three touching isotropic ellipsoids of
! materials 1,2,3
! = 'GAUSS_SPH' gaussian sphere target
! = 'HEXGONPBC' 1-d or 2-d array of hexagonal prisms
! = 'HEX_PRISM' homogeneous hexagonal prism
! = 'LYRD_SLAB' layered slab target, with up to 4 separate
! material layers
! = 'LYRSLBPBC' 1-d or 2-d array of layered rect. slab targets,
! with up to 4 material layers
! = 'MLTBLOCKS' collection of cubic blocks defined by
! data in file 'blocks.par'
! = 'RCTGLBLK3' isolated target: 3 rectangular blocks with
! centers on x-axis
! = 'RCTGLPRSM' homogeneous rectangular prism
! = 'RCTGL_PBC' 1-d or 2-d array of rectangular targets
! = 'SLAB_HOLE' rectangular block with cylindrical hole
! = 'SLBHOLPBC' 1-d or 2-d array of rectangular blocks
! with cylindrical hole
! = 'SPHERES_N' multisphere target = union of N spheres
! = 'SPHRN_PBC' 1-d or 2-d array of multisphere target
! = 'SPHROID_2' two touching spheroids with symmetry axes at
! specified angle!
! = 'SPH_ANI_N' multisphere target, with spheres that
! can have different, anisotropic, composition
! = 'TETRAHDRN' regular tetrahedron
! = 'TRILYRPBC' periodic target: 3 layer rectangular structure
! = 'TRNGLPRSM' triangular prism (homogeneous and isotropic)
! = 'UNIAXICYL' cylinder of unixaxial material
CWHERE='error reading CSHAPE from ddscat.par'
READ(IOPAR,FMT=*,ERR=99)CSHAPE
IF(CSHAPE=='aniellips')CSHAPE='ANIELLIPS'
IF(CSHAPE=='ani_ell_2')CSHAPE='ANI_ELL_2'
IF(CSHAPE=='anifilpbc')CSHAPE='ANIFILPBC'
IF(CSHAPE=='anifrmfil')CSHAPE='ANIFRMFIL'
IF(CSHAPE=='anirctngl')CSHAPE='ANIRCTNGL'
IF(CSHAPE=='bislinpbc')CSHAPE='BISLINPBC'
IF(CSHAPE=='conellips')CSHAPE='CONELLIPS'
IF(CSHAPE=='cylinder1')CSHAPE='CYLINDER1'
IF(CSHAPE=='cylndrcap')CSHAPE='CYLNDRCAP'
IF(CSHAPE=='cylndrpbc')CSHAPE='CYLNDRPBC'
IF(CSHAPE=='dskblypbc')CSHAPE='DSKBLYPBC'
IF(CSHAPE=='dskrctpbc')CSHAPE='DSKRCTPBC'
IF(CSHAPE=='dw1996tar')CSHAPE='DW1996TAR'
IF(CSHAPE=='el_in_rct')CSHAPE='EL_IN_RCT'
IF(CSHAPE=='ellipsoid')CSHAPE='ELLIPSOID'
IF(CSHAPE=='ellipso_2')CSHAPE='ELLIPSO_2'
IF(CSHAPE=='ellipso_3')CSHAPE='ELLIPSO_3'
IF(CSHAPE=='ellipspbc')CSHAPE='ELLIPSPBC'
IF(CSHAPE=='frmfilpbc')CSHAPE='FRMFILPBC'
IF(CSHAPE=='from_file')CSHAPE='FROM_FILE'
IF(CSHAPE=='gauss_sph')CSHAPE='GAUSS_SPH'
IF(CSHAPE=='hexgonpbc')CSHAPE='HEXGONPBC'
IF(CSHAPE=='hex_prism')CSHAPE='HEX_PRISM'
IF(CSHAPE=='layrdslab')CSHAPE='LAYRDSLAB'
IF(CSHAPE=='lyrslbpbc')CSHAPE='LYRSLBPBC'
IF(CSHAPE=='mltblocks')CSHAPE='MLTBLOCKS'
IF(CSHAPE=='rctglblk3')CSHAPE='RCTGLBLK3'
IF(CSHAPE=='rctglprsm')CSHAPE='RCTGLPRSM'
IF(CSHAPE=='rctgl_pbc')CSHAPE='RCTGL_PBC'
IF(CSHAPE=='rctg_rctg')CSHAPE='RCTG_RCTG'
IF(CSHAPE=='recrecpbc')CSHAPE='RECRECPBC'
IF(CSHAPE=='slab_hole')CSHAPE='SLAB_HOLE'
IF(CSHAPE=='slbholpbc')CSHAPE='SLBHOLPBC'
IF(CSHAPE=='spheres_n')CSHAPE='SPHERES_N'
IF(CSHAPE=='sphrn_pbc')CSHAPE='SPHRN_PBC'
IF(CSHAPE=='sphroid_2')CSHAPE='SPHROID_2'
IF(CSHAPE=='sph_ani_n')CSHAPE='SPH_ANI_N'
IF(CSHAPE=='tetrahdrn')CSHAPE='TETRAHDRN'
IF(CSHAPE=='trilyrpbc')CSHAPE='TRILYRPBC'
IF(CSHAPE=='trnglprsm')CSHAPE='TRNGLPRSM'
IF(CSHAPE=='uniaxicyl')CSHAPE='UNIAXICYL'
IF(CSHAPE=='ANIELLIPS'.OR. &
CSHAPE=='ANI_ELL_2'.OR. &
CSHAPE=='ANI_ELL_3'.OR. &
CSHAPE=='ANIFRMFIL'.OR. &
CSHAPE=='ANIRCTNGL'.OR. &
CSHAPE=='BISLINPBC'.OR. &
CSHAPE=='CONELLIPS'.OR. &
CSHAPE=='CYLINDER1'.OR. &
CSHAPE=='CYLNDRCAP'.OR. &
CSHAPE=='CYLNDRPBC'.OR. &
CSHAPE=='DSKBLYPBC'.OR. &
CSHAPE=='DSKRCTNGL'.OR. &
CSHAPE=='DSKRCTPBC'.OR. &
CSHAPE=='DW1996TAR'.OR. &
CSHAPE=='EL_IN_RCT'.OR. &
CSHAPE=='ELLIPSOID'.OR. &
CSHAPE=='ELLIPSO_2'.OR. &
CSHAPE=='ELLIPSO_3'.OR. &
CSHAPE=='ELLIPSPBC'.OR. &
CSHAPE=='FRMFILPBC'.OR. &
CSHAPE=='FROM_FILE'.OR. &
CSHAPE=='GAUSS_SPH'.OR. &
CSHAPE=='HEXGONPBC'.OR. &
CSHAPE=='HEX_PRISM'.OR. &
CSHAPE=='LAYRDSLAB'.OR. &
CSHAPE=='LYRSLBPBC'.OR. &
CSHAPE=='MLTBLOCKS'.OR. &
CSHAPE=='RCTGLBLK3'.OR. &
CSHAPE=='RCTGLPRSM'.OR. &
CSHAPE=='RCTGL_PBC'.OR. &
CSHAPE=='RCTG_RCTG'.OR. &
CSHAPE=='RECRECPBC'.OR. &
CSHAPE=='SLAB_HOLE'.OR. &
CSHAPE=='SLBHOLPBC'.OR. &
CSHAPE=='SPHERES_N'.OR. &
CSHAPE=='SPHRN_PBC'.OR. &
CSHAPE=='SPHROID_2'.OR. &
CSHAPE=='SPH_ANI_N'.OR. &
CSHAPE=='TETRAHDRN'.OR. &
CSHAPE=='TRILYRPBC'.OR. &
CSHAPE=='TRNGLPRSM'.OR. &
CSHAPE=='UNIAXICYL')THEN
WRITE(CMSGNM,FMT='(A,A)')CSHAPE,' - Shape definition '
CALL WRIMSG('REAPAR',CMSGNM)
ELSE
WRITE(CMSGNM,FMT='(A,A)')'CSHAPE=',CSHAPE
CALL WRIMSG('REAPAR',CMSGNM)
CALL ERRMSG('FATAL','REAPAR', &
' Unrecognized shape directive')
ENDIF
! Read shape parameters
! targets needing 0 shape parameters (skip one line):
IF(CSHAPE=='FROM_FILE'.OR.CSHAPE=='ANIFRMFIL'.OR. &
CSHAPE=='MLTBLOCKS')THEN
READ(IOPAR,FMT=*,ERR=99)
! targets needing 1 shape parameter:
ELSEIF(CSHAPE=='DW1996TAR'.OR.CSHAPE=='TETRAHDRN')THEN
CWHERE='error reading 1 shape parameter in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)SHPAR(1)
WRITE(CMSGNM,FMT='(0PF10.5,A)')SHPAR(1),' = SHPAR(1)'
CALL WRIMSG('REAPAR',CMSGNM)
! targets needing 2 shape parameters:
ELSEIF(CSHAPE=='CYLNDRCAP'.OR.CSHAPE=='UNIAXICYL')THEN
CWHERE='error reading 2 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)SHPAR(1),SHPAR(2)
DO J=1,2
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 2 shape parameters and file name:
ELSEIF(CSHAPE=='FRMFILPBC'.OR.CSHAPE=='SPHERES_N'.OR. &
CSHAPE=='SPH_ANI_N')THEN
CWHERE='error reading 2 shape parameters and filename in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)SHPAR(1),SHPAR(2),CFLSHP
CALL WRIMSG('REAPAR','shape file=')
CALL WRIMSG('REAPAR',CFLSHP)
DO J=1,2
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 3 shape parameters:
ELSEIF(CSHAPE=='ANIELLIPS'.OR.CSHAPE=='ANI_ELL_2'.OR. &
CSHAPE=='ANI_ELL_3'.OR.CSHAPE=='ANIRCTNGL'.OR. &
CSHAPE=='CYLINDER1'.OR.CSHAPE=='ELLIPSOID'.OR. &
CSHAPE=='ELLIPSO_2'.OR.CSHAPE=='ELLIPSO_3'.OR. &
CSHAPE=='HEX_PRISM'.OR.CSHAPE=='RCTGLPRSM')THEN
CWHERE='error reading 2 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)SHPAR(1),SHPAR(2),SHPAR(3)
DO J=1,3
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 3 shape parameters and file name:
ELSEIF(CSHAPE=='SPHRN_PBC')THEN
CWHERE='error reading 3 shape parameters and filename in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)SHPAR(1),SHPAR(2),SHPAR(3),CFLSHP
CALL WRIMSG('REAPAR','shape file=')
CALL WRIMSG('REAPAR',CFLSHP)
DO J=1,3
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 4 shape parameters
ELSEIF(CSHAPE=='SLAB_HOLE'.OR.CSHAPE=='TRNGLPRSM')THEN
CWHERE='error reading 4 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)SHPAR(1),SHPAR(2),SHPAR(3),SHPAR(4)
DO J=1,4
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 5 shape parameters
ELSEIF(CSHAPE=='CYLNDRPBC'.OR.CSHAPE=='DSKRCTNGL'.OR. &
CSHAPE=='HEXGONPBC'.OR.CSHAPE=='RCTGL_PBC')THEN
CWHERE='error reading 5 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)(SHPAR(J),J=1,5)
DO J=1,5
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 6 shape parameters:
ELSEIF(CSHAPE=='BISLINPBC'.OR.CSHAPE=='CONELLIPS'.OR. &
CSHAPE=='EL_IN_RCT'.OR.CSHAPE=='GAUSS_SPH'.OR. &
CSHAPE=='RCTG_RCTG'.OR.CSHAPE=='SLBHOLPBC'.OR. &
CSHAPE=='SPHROID_2')THEN
CWHERE='error reading 6 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)(SHPAR(J),J=1,6)
DO J=1,6
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 7 shape parameters:
ELSEIF(CSHAPE=='DSKRCTPBC'.OR.CSHAPE=='LAYRDSLAB')THEN
CWHERE='error reading 7 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)(SHPAR(J),J=1,7)
DO J=1,7
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 8 shape parameters:
ELSEIF(CSHAPE=='DSKBLYPBC'.OR.CSHAPE=='RECRECPBC')THEN
CWHERE='error reading 8 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)(SHPAR(J),J=1,8)
DO J=1,8
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 9 shape parameters:
ELSEIF(CSHAPE=='LYRSLBPBC'.OR.CSHAPE=='RCTGLBLK3')THEN
CWHERE='error reading 9 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)(SHPAR(J),J=1,9)
DO J=1,9
WRITE(CMSGNM,FMT='(0PF10.5,A,I1,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
! targets needing 10 shape parameters: none
! targets needing 11 shape parameters:
ELSEIF(CSHAPE=='TRILYRPBC')THEN
CWHERE='error reading 11 shape parameters in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)(SHPAR(J),J=1,11)
DO J=1,11
WRITE(CMSGNM,FMT='(0PF10.5,A,I2,A)')SHPAR(J),' = SHPAR(',J,')'
CALL WRIMSG('REAPAR',CMSGNM)
ENDDO
ELSE
CALL ERRMSG('FATAL','REAPAR', &
' No instructions for reading shape parameters for this shape')
ENDIF
! JPBC=0 if PBC are not used
! 1 if PBC in y direction only
! 2 if PBC in z direction only
! 3 if PBC in y and z directions
JPBC=0
! need to set JPBC for PBC targets
IF(CSHAPE=='BISLINPBC')THEN
IF(SHPAR(6)>0._WP)JPBC=1
JPBC=JPBC+2
ENDIF
IF(CSHAPE=='ANIFILPBC')THEN
IF(SHPAR(1)>0._WP)JPBC=1
IF(SHPAR(2)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='CYLNDRPBC')THEN
IF(SHPAR(4)>0._WP)JPBC=1
IF(SHPAR(5)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='DSKBLYPBC')THEN
IF(SHPAR(7)>0._WP)JPBC=1
IF(SHPAR(8)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='DSKRCTPBC')THEN
IF(SHPAR(6)>0._WP)JPBC=1
IF(SHPAR(7)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='FRMFILPBC')THEN
IF(SHPAR(1)>0._WP)JPBC=1
IF(SHPAR(2)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='HEXGONPBC')THEN
IF(SHPAR(4)>0._WP)JPBC=1
IF(SHPAR(5)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='LYRSLBPBC')THEN
IF(SHPAR(8)>0._WP)JPBC=1
IF(SHPAR(9)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='RCTGL_PBC')THEN
IF(SHPAR(4)>0._WP)JPBC=1
IF(SHPAR(5)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='RECRECPBC')THEN
IF(SHPAR(7)>0._WP)JPBC=1
IF(SHPAR(8)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='SLBHOLPBC')THEN
IF(SHPAR(5)>0._WP)JPBC=1
IF(SHPAR(6)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='SPHRN_PBC')THEN
IF(SHPAR(2)>0._WP)JPBC=1
IF(SHPAR(3)>0._WP)JPBC=JPBC+2
ENDIF
IF(CSHAPE=='TRILYRPBC')THEN
IF(SHPAR(10)>0._WP)JPBC=1
IF(SHPAR(11)>0._WP)JPBC=JPBC+2
ENDIF
!=======================================================================
! following code disabled 03.01.29
! (retain for use in future noncubic version)
! Obtain lattice anisotropy parameters DX(1-3)
! For cubic lattice, DX(1)=DX(2)=DX(3)
! Note that we do not require here that DX(1)*DX(2)*DX(3)=1 upon
! input; DX is renormalized here before being returned to DDSCAT
! READ(IOPAR,FMT=*,ERR=99)DX(1),DX(2),DX(3)
! WRITE(CMSGNM,FMT='(F8.3,F8.3,F8.3,A)')DX,
! & ' = relative lattice spacings dx,dy,dz'
! CALL WRIMSG('REAPAR',CMSGNM)
! DELTA=(DX(1)*DX(2)*DX(3))**(1./3.)
! DX(1)=DX(1)/DELTA
! DX(2)=DX(2)/DELTA
! DX(3)=DX(3)/DELTA
! WRITE(CMSGNM,FMT='(F8.3,F8.3,F8.3,A)')DX,
! & ' = normalized lattice spacings dx,dy,dz'
! CALL WRIMSG('REAPAR',CMSGNM)
! and replaced by following:
DX(1)=1._WP
DX(2)=1._WP
DX(3)=1._WP
!=======================================================================
! Obtain names of file(s) containing dielectric function(s)
! NCOMP = number of different dielectric functions
! CFLEPS(1-NCOMP) = names of files containing dielectric functions
CWHERE='error reading NCOMP in ddscat.par'
READ(IOPAR,FMT=*,ERR=99)NCOMP
WRITE(CMSGNM,FMT='(A,I4)')' NCOMP=',NCOMP
CALL WRIMSG('REAPAR',CMSGNM)
! Check that MXCOMP is large enough
IF(NCOMP>MXCOMP)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP > MXCOMP: must increase MXCOMP in DDSCAT')
!*** Check that NCOMP=2 if CSHAPE=UNIAXICYL *******************************
! 3 ANIELLIPS
! 2 ELLIPSO_2
! 3 ELLIPSO_3
! 6 ANI_ELL_2
! 9 ANI_ELL_3
! 2 CONELLIPS
! 2 SPHROID_2
! 2-4 LYRD_SLAB
! 1-4 LYRSLBPBC
IF(CSHAPE=='UNIAXICYL'.AND.NCOMP/=2)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP must be 2 for option UNIAXICYL')
IF(CSHAPE=='ANIELLIPS'.AND.NCOMP/=3)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP must be 3 for option ANIELLIPS')
IF(CSHAPE=='ELLIPSO_2'.AND.NCOMP/=2)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP must be 2 for option ELLIPSO_2')
IF(CSHAPE=='ELLIPSO_3'.AND.NCOMP/=3)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP must be 3 for option ELLIPSO_3')
IF(CSHAPE=='ANI_ELL_2'.AND.NCOMP/=6)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP must be 6 for option ANI_ELL_2')
IF(CSHAPE=='ANI_ELL_3'.AND.NCOMP/=9)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP must be 9 for option ANI_ELL_3')
IF(CSHAPE=='CONELLIPS'.AND.NCOMP/=2)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP must be 2 for option CONELLIPS')
IF(CSHAPE=='SPHROID_2'.AND.NCOMP/=2)CALL ERRMSG('FATAL','REAPAR', &
'NCOMP must be 2 for option SPHROID_2')
IF(CSHAPE=='LYRSLBPBC')THEN
J=4
IF(SHPAR(7)<=0._WP)J=3
IF(SHPAR(7)<=0._WP.AND.SHPAR(6)<=0._WP)J=2