      PROGRAM CYLABS
C***********************************************************************
C*                                                                     *
C* CYLABS: a FORTRAN77 program for calculating X-ray absorption for    *
C* a cylindric single crystal sample.                                  *
C*                                                                     *
C*                           by M.NARDELLI                             *
C*                                                                     *
C* Istituto di Chimica Generale ed Inorganica della Universita` degli  *
C* Studi di Parma, Centro di Studio per la Strutturistica Diffrattome- *
C* trica del C.N.R., Viale delle Scienze, I-43100 Parma, Italy.        *
C*                                                                     *
C* The program calculates: the absorption correction by interpolating  *
C* data from a table giving the absorption correction factors. The     *
C* table is that for cylindric crystals: Int.Tab. vol. C (1992) p. 523;*
C* Tab. 6.3.3.3 for values of muR in the range 0-2.5, and Int.Tab.     *
C* vol.II (1959) p.302; Table 5.3.6 B for muR in the range 2.6-10.0.   *
C* The interpolation procedure illustrated by C.W.Dwiggins jr (Acta    *
C* Cryst.(1975).A31,146-148) is used with some modification.           *
C* The program corrects for absorption the observed squared structure  *
C* factors or structure factors in the input table for SHELXL-93 or    *
C* SHELX- 76, respectively (FORMAT(3I4,2F8.2,3X,I1).                   *
C* The muR value, the absorption correction as a function of theta,    *
C* the max., min. and mean values of transmission factors, error       *
C* messages (if any) are quoted in a file named CYLABS.RIS             *
C* The computer interactively asks :                                   *
C*         the name of the input file (the same as for SHELX76)        *
C*         the name of the output file (input for SHELX76)             *
C*         the name of the compound                                    *
C*         a,b,c,alpha,beta,gamma                                      *
C*         the muR value                                               *
C*         the code for F^2 or F                                       *
C*         the value of lambda (wave-length of the radiation)          *
C* Subroutines used:                                                   *
C*         RCCNST                                                      *
C*         COSTK(X1,X2,X3,Y1,Y2,Y3,AK0,AK1,AK2)                        *
C*         COSTL(T1,T2,T3,T4,Y1,Y2,Y3,Y4,AL0,AL1,AL2,AL3)              *
C***********************************************************************
      CHARACTER*1 TIT(80)
      CHARACTER*14 FILIN,FILOUT
      INTEGER H
      DOUBLE PRECISION ST,AL0,AL1,AL2,AL3,SN(3),CS(3),CR(3),CSR(3)
      DIMENSION C(3),A(3),AF(101,19),RM(101),
     1 TH(19),AS(19),AK(3,19)
      COMMON/GRP1/RD
     1      /GRP2/C,A
     2      /GRP3/SN,CS,CR,CSR
      DATA((AF(I,K),I=1,101),K=1,2)/1.0,1.1843,1.4009,1.6548,1.9522,
     1 2.2996,2.7047,3.1762,3.7236,4.3578,5.0907,5.9356,6.907,8.021,
     2 9.294,10.746,12.397,14.267,16.379,18.76,21.43,24.41,27.74,
     3 31.44,35.54,40.06,44.7,50.1,56.1,62.5,69.5,77.2,85.4,94.2,
     4 104,114,125,136,149,162,175,190,206,222,239,257,275,295,316,
     5 337,359,383,407,432,458,485,513,542,573,604,636,670,704,740,
     6 777,814,853,894,935,978,1022,1067,1113,1161,1210,1260,1311,
     7 1364,1418,1474,1530,1597,1664,1731,1798,1865,1932,1999,2066,
     8 2133,2200,2285,2370,2455,2540,2625,2710,2795,2880,2965,3050,
     9 1.0,1.1843,1.4007,1.6544,1.9513,2.2979,2.7017,3.1712,3.7157,
     1 4.3456,5.0724,5.9089,6.869,7.967,9.219,10.643,12.257,14.080,
     2 16.131,18.43,21.00,23.87,27.04,30.55,34.41,38.65,43.1,48.1,
     3 53.5,59.4,65.8,72.6,79.9,87.6,95.9,105,114,124,134,145,156,
     4 168,180,193,206,220,234,249,264,280,296,313,330,348,366,384,
     5 403,422,442,462,483,504,525,547,569,591,614,638,661,685,710,
     6 735,760,785,811,837,864,891,918,945,974,1004,1033,1063,1092,
     7 1122,1152,1181,1211,1240,1270,1303,1336,1369,1402,1435,1468,
     8 1501,1534,1567,1600/
      DATA((AF(I,K),I=1,101),K=3,4)/1.0,1.1842,1.4002,1.6531,1.9485,
     1 2.2926,2.6926,3.1561,3.6919,4.3093,5.0185,5.8305,6.757,7.810,
     2 9.003,10.349,11.862,13.555,15.441,17.53,19.84,22.39,25.17,28.20,
     3 31.49,35.05,38.7,42.8,47.2,51.8,56.7,61.8,67.3,72.9,78.9,85.0,
     4 91.4,98.0,105,112,119,127,134,142,150,158,166,175,183,192,200,
     5 209,218,228,237,246,255,265,275,284,294,304,314,324,334,345,355,
     6 365,376,386,397,408,418,429,440,451,462,473,484,495,507,519,530,
     7 542,553,565,577,588,600,611,623,635,647,659,671,683,695,707,719,
     8 731,743,1.0,1.1840,1.3995,1.6510,1.9439,2.2840,2.6775,3.1315,
     9 3.6532,4.2507,4.9323,5.7065,6.582,7.568,8.674,9.907,11.276,
     1 12.788,14.450,16.267,18.24,20.38,22.69,25.16,27.79,30.59,33.4,
     2 36.5,39.7,43.1,46.6,50.3,54.0,57.9,61.9,65.9,70.1,74.4,78.7,83.1,
     3 87.6,92.1,96.7,101,106,111,116,121,125,130,135,140,145,150,156,
     4 161,166,171,176,182,187,192,198,203,209,214,219,225,231,236,242,
     5 247,253,259,264,270,276,281,287,293,299,305,311,316,322,328,334,
     6 340,345,351,357,363,369,375,381,387,393,399,405,411,417/
      DATA((AF(I,K),I=1,101),K=5,6)/1.0,1.1838,1.3984,1.6481,1.9376,
     1 2.2721,2.6570,3.0982,3.6015,4.1733,4.8196,5.5466,6.360,7.266,
     2 8.268,9.372,10.581,11.897,13.323,14.858,16.50,18.25,20.11,22.07,
     3 24.13,26.28,28.5,30.8,33.2,35.7,38.2,40.8,43.5,46.2,48.9,51.7,
     4 54.6,57.4,60.4,63.3,66.3,69.3,72.3,75.4,78.5,81.6,84.7,87.8,
     5 91.0,94.2,97.4,101,104,107,110,114,117,120,124,127,130,134,137,
     6 140,144,147,151,154,158,161,165,168,172,175,179,182,186,190,193,
     7 197,200,204,207,211,214,218,222,225,229,232,236,240,243,247,250,
     8 254,258,261,265,268,272,1.0,1.1835,1.3970,1.6443,1.9296,2.2572,
     9 2.6317,3.0575,3.5392,4.0812,4.6877,5.3624,6.109,6.929,7.826,
     1 8.800,9.852,10.982,12.189,13.470,14.824,16.247,17.74,19.29,20.90,
     2 22.56,24.2,26.0,27.8,29.6,31.5,33.4,35.4,37.3,39.3,41.3,43.3,
     3 45.4,47.5,49.5,51.6,53.8,55.9,58.0,60.2,62.3,64.5,66.7,68.9,71.1,
     4 73.3,75.5,77.8,80.0,82.3,84.5,86.8,89.1,91.4,93.6,95.9,98.3,101,
     5 103,105,108,110,112,115,117,119,122,124,127,129,131,134,136,139,
     6 141,144,146,149,151,153,156,158,160,162,165,167,169,172,174,177,
     7 179,181,184,186,189,191/
      DATA((AF(I,K),I=1,101),K=7,8)/1.0,1.1832,1.3953,1.6398,1.9201,
     1 2.2398,2.6023,3.0111,3.4691,3.9792,4.5439,5.1649,5.8436,6.581,
     2 7.376,8.320,9.141,10.106,11.125,12.194,13.311,14.472,15.675,
     3 16.92,18.19,19.50,20.8,22.2,23.6,25.0,26.4,27.9,29.3,30.8,32.3,
     4 33.8,35.3,36.9,38.4,39.9,41.5,43.1,44.6,46.2,47.8,49.4,51.0,52.6,
     5 54.2,55.8,57.5,59.1,60.7,62.4,64.0,65.7,67.3,69.0,70.7,72.4,74.0,
     6 75.7,77.4,79.1,80.8,82.5,84.2,85.9,87.6,89.4,91.1,92.8,94.6,96.3,
     7 98.1,99.8,102,103,105,107,109,111,112,114,116,118,119,121,123,
     8 124,126,128,129,131,133,134,136,138,140,141,143,1.0,1.1828,
     9 1.3934,1.6347,1.9094,2.2204,2.5701,2.9607,3.3941,3.8718,4.3948,
     1 4.9636,5.5782,6.238,6.942,7.689,8.477,9.304,10.168,11.066,11.995,
     2 12.953,13.938,14.947,15.978,17.03,18.1,19.2,20.3,21.4,22.5,23.7,
     3 24.8,26.0,27.1,28.3,29.5,30.7,31.9,33.1,34.3,35.5,36.7,37.9,39.1,
     4 40.4,41.6,42.8,44.1,45.3,46.6,47.8,49.1,50.3,51.6,52.9,54.1,55.4,
     5 56.7,58.0,59.2,60.5,61.8,63.1,64.4,65.7,67.0,68.3,69.6,70.9,72.2,
     6 73.6,74.9,76.2,77.5,78.9,80.2,81.5,82.9,84.2,85.5,86.8,88.0,89.3,
     7 90.5,91.8,93.1,94.3,95.6,96.8,98.1,99.4,101,102,103,105,106,107,
     8 108,110,111/
      DATA((AF(I,K),I=1,101),K=9,10)/1.0,1.1823,1.3912,1.6290,1.8979,
     1 2.1996,2.5359,2.9081,3.3169,3.7629,4.4261,4.7660,5.3219,5.9125,
     2 6.536,7.192,7.877,8.589,9.327,10.089,10.871,11.673,12.493,13.328,
     3 14.177,15.040,15.9,16.8,17.7,18.6,19.5,20.4,21.4,22.3,23.2,24.2,
     4 25.1,26.1,27.0,28.0,28.9,29.9,30.9,31.8,32.8,33.8,34.8,35.8,36.7,
     5 37.7,38.7,39.7,40.7,41.7,42.7,43.7,44.7,45.7,46.8,47.8,48.8,49.8,
     6 50.8,51.9,52.9,53.9,54.9,56.0,57.0,58.0,59.1,60.1,61.2,62.2,63.3,
     7 64.3,65.4,66.4,67.5,68.5,69.6,70.6,71.6,72.6,73.6,74.6,75.6,76.6,
     8 77.6,78.6,79.6,80.6,81.6,82.6,83.6,84.7,85.7,86.7,87.7,88.7,89.7,
     7 1.0,1.1818,1.3889,1.6230,1.8857,2.1781,2.5010,2.8549,3.2400,
     8 3.6560,4.1022,4.5776,5.0811,5.6110,6.166,6.744,7.344,7.963,8.600,
     9 9.253,9.921,10.602,11.295,11.999,12.711,13.433,14.2,14.9,15.6,
     1 16.4,17.1,17.9,18.7,19.4,20.2,21.0,21.7,22.5,23.3,24.1,24.9,25.7,
     2 26.5,27.3,28.1,28.9,29.7,30.5,31.3,32.1,32.9,33.7,34.5,35.3,36.1,
     3 37.0,37.8,38.6,39.4,40.3,41.1,41.9,42.7,43.6,44.4,45.2,46.1,46.9,
     4 47.8,48.6,49.4,50.3,51.1,52.0,52.8,53.7,54.5,55.4,56.2,57.1,58.0,
     5 58.8,59.5,60.3,61.0,61.8,62.6,63.3,64.1,64.8,65.6,66.5,67.3,68.2,
     6 69.0,69.9,70.8,71.6,72.5,73.3,74.2/
      DATA((AF(I,K),I=1,101),K=11,12)/1.0,1.1813,1.3865,1.6169,1.8733,
     1 2.1564,2.4662,2.8028,3.1656,3.5538,3.9664,4.4022,4.8598,5.3376,
     2 5.8341,6.348,6.877,7.420,7.976,8.544,9.122,9.709,10.304,10.906,
     3 11.515,12.130,12.8,13.4,14.0,14.6,15.3,15.9,16.6,17.2,17.9,18.5,
     4 19.2,19.8,20.5,21.2,21.8,22.5,23.2,23.8,24.5,25.2,25.8,26.5,27.2,
     5 27.9,28.6,29.2,29.9,30.6,31.3,32.0,32.7,33.4,34.0,34.7,35.4,36.1,
     6 36.8,37.5,38.2,38.9,39.6,40.3,41.0,41.7,42.4,43.1,43.8,44.6,45.3,
     7 46.0,46.7,47.4,48.1,48.9,49.5,50.2,50.9,51.5,52.2,52.9,53.6,54.3,
     8 54.9,55.6,56.3,57.0,57.6,58.3,59.0,59.7,60.3,61.0,61.7,62.3,63.0,
     9 1.0,1.1808,1.3841,1.6108,1.8611,2.1352,2.4327,2.7530,3.0953,
     1 3.4584,3.8413,4.2424,4.6604,5.0938,5.5413,6.002,6.473,6.955,
     2 7.446,7.946,8.452,8.965,9.484,10.008,10.537,11.069,11.6,12.2,
     3 12.7,13.2,13.8,14.3,14.9,15.5,16.0,16.6,17.1,17.7,18.3,18.8,19.4,
     4 20.0,20.6,21.1,21.7,22.3,22.9,23.4,24.0,24.6,25.2,25.8,26.3,26.9,
     5 27.5,28.1,28.7,29.3,29.9,30.5,31.0,31.6,32.2,32.8,33.4,34.0,34.6,
     6 35.2,35.8,36.4,37.0,37.6,38.2,38.8,39.5,40.1,40.7,41.3,41.9,42.5,
     7 43.1,43.7,44.2,44.8,45.3,45.9,46.5,47.0,47.6,48.1,48.7,49.3,49.8,
     8 50.4,50.9,51.5,52.1,52.6,53.2,53.7,54.3/
      DATA((AF(I,K),I=1,101),K=13,14)/1.0,1.1802,1.3818,1.6049,1.8495,
     1 2.1152,2.4012,2.7068,3.0307,3.3717,3.7286,4.0998,4.4842,4.8805,
     2 5.2873,5.7036,6.128,6.561,7.000,7.444,7.895,8.349,8.808,9.271,
     3 9.736,10.205,10.7,11.2,11.7,12.1,12.6,13.1,13.6,14.1,14.6,15.1,
     4 15.6,16.1,16.6,17.1,17.6,18.1,18.6,19.1,19.6,20.1,20.6,21.1,21.6,
     5 22.1,22.6,23.1,23.6,24.1,24.6,25.2,25.7,26.2,26.7,27.2,27.7,28.3,
     6 28.8,29.3,29.8,30.3,30.9,31.4,31.9,32.4,33.0,33.5,34.0,34.6,35.1,
     7 35.6,36.1,36.7,37.2,37.8,38.3,38.8,39.3,39.8,40.3,40.7,41.2,41.7,
     8 42.2,42.7,43.2,43.7,44.2,44.7,45.2,45.7,46.2,46.7,47.2,47.7,48.2,
     9 1.0,1.1798,1.3796,1.5994,1.8388,2.0969,2.3728,2.6653,2.9732,
     1 3.2951,3.6298,3.9759,4.3322,4.6976,5.0711,5.4516,5.8385,6.231,
     2 6.628,7.030,7.435,7.843,8.255,8.669,9.086,9.505,9.93,10.4,10.8,
     3 11.2,11.6,12.1,12.5,12.9,13.4,13.8,14.2,14.7,15.1,15.6,16.0,16.4,
     4 16.9,17.3,17.8,18.2,18.7,19.1,19.6,20.0,20.5,20.9,21.4,21.8,22.3,
     5 22.7,23.2,23.6,24.1,24.5,25.0,25.4,25.9,26.4,26.8,27.3,27.7,28.2,
     6 28.7,29.1,29.6,30.1,30.5,31.0,31.5,31.9,32.4,32.9,33.3,33.8,34.3,
     7 34.7,35.2,35.6,36.0,36.4,36.9,37.3,37.7,38.2,38.6,39.0,39.5,39.9,
     8 40.4,40.8,41.2,41.7,42.1,42.6,43.0/
      DATA((AF(I,K),I=1,101),K=15,16)/1.0,1.1793,1.3777,1.5946,1.8293,
     1 2.0809,2.3480,2.6295,2.9239,3.2299,3.5462,3.8717,4.2051,4.5456,
     2 4.8922,5.2441,5.6007,5.961,6.326,6.693,7.063,7.436,7.810,8.187,
     3 8.565,8.945,9.33,9.72,10.1,10.5,10.9,11.3,11.7,12.0,12.4,12.8,
     4 13.2,13.6,14.0,14.4,14.8,15.2,15.6,16.0,16.4,16.8,17.2,17.6,18.0,
     5 18.4,18.8,19.2,19.6,20.0,20.4,20.9,21.3,21.7,22.1,22.5,22.9,23.3,
     6 23.7,24.1,24.6,25.0,25.4,25.8,26.2,26.6,27.0,27.5,27.9,28.3,28.7,
     7 29.2,29.6,30.0,30.4,30.8,31.3,31.7,32.1,32.5,32.9,33.2,33.6,34.0,
     8 34.4,34.8,35.2,35.6,36.0,36.4,36.8,37.2,37.5,37.9,38.3,38.7,39.1,
     9 1.0,1.1790,1.3761,1.5906,1.8215,2.0677,2.3277,2.6003,2.8839,
     1 3.1772,3.4790,3.7882,4.1038,4.4248,4.7506,5.0804,5.4136,5.4799,
     2 6.089,6.430,6.773,7.118,7.464,7.812,8.161,8.511,8.88,9.23,9.59,
     3 9.95,10.3,10.7,11.0,11.4,11.7,12.1,12.5,12.8,13.2,13.6,13.9,
     3 14.3,14.6,15.0,15.4,15.7,16.1,16.5,16.9,17.2,17.6,18.0,18.3,18.7,
     4 19.1,19.5,19.8,20.2,20.6,21.0,21.3,21.7,22.1,22.5,22.8,23.2,23.6,
     5 24.0,24.4,24.8,25.1,25.5,25.9,26.3,26.7,27.1,27.4,27.8,28.2,28.6,
     6 29.0,29.4,29.7,30.0,30.4,30.7,31.1,31.5,31.8,32.2,32.5,32.9,33.2,
     7 33.6,33.9,34.3,34.7,35.0,35.4,35.7,36.1/
      DATA((AF(I,K),I=1,101),K=17,19)/1.0,1.1787,1.3749,1.5876,1.8157,
     1 2.0579,2.3126,2.5786,2.8542,3.1383,3.4295,3.7269,4.0295,4.3365,
     2 4.6471,4.6909,5.2773,5.5960,5.9166,6.239,6.562,6.887,7.213,7.540,
     3 7.868,8.196,8.53,8.87,9.20,9.53,9.86,10.2,10.5,10.9,11.2,11.5,
     4 11.9,12.2,12.6,12.9,13.2,13.6,13.9,14.3,14.6,14.9,15.3,15.6,16.0,
     5 16.3,16.7,17.0,17.3,17.7,18.0,18.4,18.7,19.1,19.4,19.8,20.1,20.5,
     6 20.8,21.2,21.5,21.9,22.2,22.6,23.0,23.3,23.7,24.0,24.4,24.7,25.1,
     7 25.5,25.8,26.2,26.5,26.9,27.3,27.6,28.0,28.3,28.6,28.9,29.3,29.6,
     8 29.9,30.3,30.6,30.9,31.3,31.6,31.9,32.2,32.6,32.9,33.2,33.6,33.9,
     9 1.0,1.1785,1.3741,1.5857,1.8121,2.0518,2.3033,2.5651,2.8359,
     1 3.1142,3.3990,3.6891,3.9838,4.2821,4.5835,4.8875,5.1935,5.5014,
     2 5.8107,6.121,6.433,6.745,7.059,7.372,7.687,8.002,8.33,8.64,8.96,
     3 9.28,9.60,9.92,10.2,10.6,10.9,11.2,11.5,11.8,12.2,12.5,12.8,13.1,
     4 13.5,13.8,14.1,14.4,14.8,15.1,15.4,15.7,16.1,16.4,16.7,17.1,17.4,
     5 17.7,18.0,18.4,18.7,19.0,19.4,19.7,20.0,20.4,20.7,21.0,21.4,21.7,
     6 22.1,22.4,22.7,23.1,23.4,23.8,24.1,24.4,24.8,25.1,25.5,25.8,26.1,
     7 26.4,26.7,27.0,27.3,27.6,28.0,28.3,28.6,28.9,29.2,29.5,29.8,30.2,
     8 30.5,30.8,31.1,31.4,31.8,32.1,32.4,1.0,1.1785,1.3739,1.5851,
     7 1.8108,2.0497,2.3001,2.5606,2.8297,3.1061,3.3886,3.6763,3.9682,
     8 4.2636,4.5619,4.8625,5.1650,5.4691,5.7746,6.081,6.389,6.697,
     9 7.006,7.315,7.625,7.935,8.27,8.58,8.90,9.21,9.53,9.84,10.2,10.5,
     1 10.8,11.1,11.4,11.7,12.1,12.4,12.7,13.0,13.3,13.7,14.0,14.3,
     2 14.6,14.9,15.3,15.6,15.9,16.2,16.6,16.9,17.2,17.5,17.9,18.2,18.5,
     3 18.8,19.2,19.5,19.8,20.2,20.5,20.8,21.1,21.5,21.8,22.1,22.5,22.8,
     4 23.1,23.5,23.8,24.1,24.5,24.8,25.2,25.5,25.8,26.1,26.4,26.7,27.0,
     5 27.3,27.7,28.0,28.3,28.6,28.9,29.2,29.5,29.8,30.1,30.4,30.7,31.0,
     6 31.3,31.6,31.9/
      DATA RM/ 0.0,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,
     1 1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.1,2.2,2.3,2.4,2.5,2.6,2.7,2.8,2.9,
     2 3.0,3.1,3.2,3.3,3.4,3.5,3.6,3.7,3.8,3.9,4.0,4.1,4.2,4.3,4.4,
     3 4.5,4.6,4.7,4.8,4.9,5.0,5.1,5.2,5.3,5.4,5.5,5.6,5.7,5.8,5.9,6.0,
     4 6.1,6.2,6.3,6.4,6.5,6.6,6.7,6.8,6.9,7.0,7.1,7.2,7.3,7.4,7.5,7.6,
     5 7.7,7.8,7.9,8.0,8.1,8.2,8.3,8.4,8.5,8.6,8.7,8.8,8.9,9.0,9.1,9.2,
     6 9.3,9.4,9.5,9.6,9.7,9.8,9.9,10.0/
      DATA TH/.0,.08727,.17453,.26180,.34907,.43633,.52360,.61087,
     1 .69813,.78540,.87266,0.95993,1.04720,1.13446,1.22173,1.30900,
     2 1.39626,1.48353,1.57080/
      RD=0.017453293
      WRITE(*,100)
100   FORMAT(' Key in the input file'/' ?>')
      READ(*,'(A14)')FILIN
      WRITE(*,101)
101   FORMAT(' Key in the output file'/' ?>')
      READ(*,'(A14)')FILOUT
      OPEN(UNIT=13,FILE=FILIN,FORM='FORMATTED')
      REWIND(UNIT=13)
      OPEN(UNIT=15,FILE=FILOUT,FORM='FORMATTED')
      REWIND(UNIT=15)
      OPEN(UNIT=16,FILE='CYLABS.RIS',FORM='FORMATTED')
      REWIND(UNIT=16)
C-----Read title, parameters
      WRITE(*,102)
102   FORMAT(' Key in the name of the compound'/' ?>')
      READ(*,'(80A1)') TIT
      WRITE(16,'(1X,80A1)') TIT
      WRITE(*,103)
103   FORMAT(' Key in a,b,c,alpha,beta,gamma'/' ?: ')
      READ(*,*) C,A
C-----Calculate reciprocal constants
      CALL RCCNST
C-----Read the muR value and choose the standard muR values for
C-----interpolation
      WRITE(*,104)
104   FORMAT(' Key in the muR value'/' ?: ')
      READ(*,*)X0
      WRITE(*,111)
111   FORMAT(' Key in the code for structure factors'/
     1 '   0  if F^2'/'   1  if F'/' ?:>')
      READ(*,*)KF
      IF(X0.LT.0.1.OR.X0.GT.10) GO TO 1001
      DO 1 I=1,101
      IF(X0.EQ.RM(I+1)) GO TO 2
      IF(X0.GT.RM(I+1)) GO TO 1
      X1=RM(I-1)
      X2=RM(I)
      X3=RM(I+1)
      DO 3 K=1,19
      Y1=AF(I-1,K)
      Y2=AF(I,K)
      Y3=AF(I+1,K)
      CALL COSTK(X1,X2,X3,Y1,Y2,Y3,AK0,AK1,AK2)
3     AS(K)=EXP(AK0+AK1*X0+AK2*X0**2)
      GO TO 9
1     CONTINUE
      WRITE(15,'(''There is something wrong with your data'')')
      STOP
2     DO 6 K=1,19
6     AS(K)=AF(I+1,K)
9     WRITE(16,'('' muR = '',F7.4/'' Theta'',9X,''A*'')')X0
      SAS=0.0
      DO 200 I=1,19
      SAS=SAS+AS(I)
      THD=TH(I)/RD
200   WRITE(16,'(1X,F5.2,5X,F9.4)')THD,AS(I)
      SA=SAS/19
C-----Read wavelength
      WRITE(*,105)
105   FORMAT(' Key in the wavelength'/' ?: ')
      READ(*,*) WL
C      IF(KF.EQ.0)WRITE(16,110)
C110   FORMAT(/18X,'Fo^2',4X,'Fo^2(corr.)'/)
C      IF(KF.EQ.1)WRITE(16,112)
C112   FORMAT(/16X,'Fo',1X,'Fo(corr.)'/)
C-----Read indices and structure factors
      TFMAX=0.0
      TFMIN=10000.0
      NN=0
      SS=0
      TH(1)=0.043633
      AS0=AS(1)
      AS(1)=(AS(1)+AS(2))/2.0
1000  CONTINUE
      READ(13,106)H,K,L,FO,SF,IN
106   FORMAT(3I4,2F8.2,3X,I1)
      IF(H.EQ.0.AND.K.EQ.0.AND.L.EQ.0) GO TO 1002
C-----Calculate sin(theta) and theta
      SNT=(WL/2.)*DSQRT((H*CR(1))**2+(K*CR(2))**2+(L*CR(3))**2+
     1        2.*(K*L*CR(2)*CR(3)*CSR(1)+
     2            L*H*CR(3)*CR(1)*CSR(2)+
     3            H*K*CR(1)*CR(2)*CSR(3)))
      IF(SNT.GT.1.0) WRITE(16,504)SNT
504   FORMAT(' sin(T0)=',F8.5)
      IF(SNT.GT.1.0) GO TO 1000
      T0=ASIN(SNT)
      ST=SNT
C-----Choose the standard theta values for interpolation
      DO 10 I=1,19
      KK=I
      IF(ABS(T0-TH(I)).LE.0.0001) GO TO 15
      IF(T0.GT.TH(2))GO TO 11
      T1=TH(1)
      T2=TH(2)
      T3=TH(3)
      T4=TH(4)
      Y1=AS(1)
      Y2=AS(2)
      Y3=AS(3)
      Y4=AS(4)
      GO TO 14
11    IF(T0.LT.TH(18))GO TO 12
      T1=TH(16)
      T2=TH(17)
      T3=TH(18)
      T4=TH(19)
      Y1=AS(16)
      Y2=AS(17)
      Y3=AS(18)
      Y4=AS(19)
      GO TO 14
12    IF(T0.LT.TH(I).OR.T0.GT.TH(I+1))GO TO 10
      IF(ABS(T0-TH(I)).LE.0.0001)THEN
      KK=I
      GO TO 15
      ENDIF
      IF(ABS(T0-TH(I+1)).LE.0.0001)THEN
      KK=I+1
      GO TO 15
      ENDIF
      IF(T0.GT.TH(2).AND.T0.LT.TH(3))THEN
      T1=TH(2)
      T2=TH(3)
      T3=TH(4)
      T4=TH(5)
      Y1=AS(2)
      Y2=AS(3)
      Y3=AS(4)
      Y4=AS(5)
      GO TO 14
      ENDIF
      T1=TH(I-1)
      T2=TH(I)
      T3=TH(I+1)
      T4=TH(I+2)
      Y1=AS(I-1)
      Y2=AS(I)
      Y3=AS(I+1)
      Y4=AS(I+2)
      GO TO 14
10    CONTINUE
14    CALL COSTL(T1,T2,T3,T4,Y1,Y2,Y3,Y4,AL0,AL1,AL2,AL3)
      ST=ST**2
      AC=(AL0+AL1*ST+AL2*ST**2+AL3*ST**3)
      GO TO 16
15    AC=AS(KK)
16    IF(KF.EQ.1)AC=SQRT(AC)
      IF(TFMAX.LT.1/AC)TFMAX=1/AC
      IF(TFMIN.GT.1/AC)TFMIN=1/AC
      T0GRD=T0/RD
      IF(AC.GT.AS0)WRITE(16,'('' AC is greater than AS(1)'')')
C     WRITE(16,'(''T0='',F6.2,3X,'' AC='',F8.4)')T0GRD,AC
      FCR=FO*AC/AS0
      SFC=SF*AC/AS0
      WRITE(15,107)H,K,L,FCR,SFC,IN
107   FORMAT(3I4,2F8.2,3X,I1)
C     WRITE(16,108)H,K,L,FO,FCR
C108  FORMAT(3I4,2F12.2)
      SS=SS+AC
      NN=NN+1
      GO TO 1000
1002  WRITE(15,106)H,K,L,FO,SF,IN
      SS=NN/SS
C     TFMAX=TFMAX/SS
C     TFMIN=TFMIN/SS
      WRITE(16,114)TFMAX,TFMIN,SS
114   FORMAT(' Transmission factor max. = ',F9.5/
     1 ' Transmission factor min. = ', F9.5/
     2 ' Transmission factor mean = ', F9.5)
      GO TO 2000
1001  WRITE(*,109)
109   FORMAT(' No calculation can be made with your muR value!')
2000  STOP
      END
 
      SUBROUTINE  RCCNST
C-----Calcola le costanti reciproche a*,b*,c*,cos(alpha*),cos(beta*),
C-----cos(gamma*) da quelle dirette.
      DOUBLE PRECISION SN(3),CS(3),CR(3),CSR(3),AD(3),S,VD,AA
      DIMENSION C(3),A(3)
      COMMON/GRP1/RD
     1      /GRP2/C,A
     2      /GRP3/SN,CS,CR,CSR
      S=0D0
      DO 1 I=1,3
      AD(I)=A(I)*RD
      SN(I)=DSIN(AD(I))
      CS(I)=DCOS(AD(I))
1     S=S+AD(I)
      S=S/2D0
      VD=2*C(1)*C(2)*C(3)
      AA=DSIN(S)
      DO 2 I=1,3
2     AA=AA*DSIN(S-AD(I))
      VD=VD*DSQRT(AA)
      CR(1)=C(2)*C(3)*SN(1)/VD
      CR(2)=C(3)*C(1)*SN(2)/VD
      CR(3)=C(1)*C(2)*SN(3)/VD
      CSR(1)=(CS(2)*CS(3)-CS(1))/(SN(2)*SN(3))
      CSR(2)=(CS(3)*CS(1)-CS(2))/(SN(3)*SN(1))
      CSR(3)=(CS(1)*CS(2)-CS(3))/(SN(1)*SN(2))
      RETURN
      END
 
      SUBROUTINE COSTK(X1,X2,X3,Y1,Y2,Y3,AK0,AK1,AK2)
C-----This routine solves a system of three linear equations
      Y1=ALOG(Y1)
      Y2=ALOG(Y2)
      Y3=ALOG(Y3)
      AA=(X2*X3**2-X3*X2**2)
      BB=(X3*X1**2-X1*X3**2)
      CC=(X1*X2**2-X2*X1**2)
      D1=AA+BB+CC
      AK0=(Y1*AA+Y2*BB+Y3*CC)/D1
      AK1=((Y1-AK0)*X2**2-(Y2-AK0)*X1**2)/CC
      AK2=(Y1-AK0-AK1*X1)/X1**2
      RETURN
      END
 
      SUBROUTINE COSTL(T1,T2,T3,T4,Y1,Y2,Y3,Y4,AL0,AL1,AL2,AL3)
C-----This routine solves a system of four linear equations
      COMMON/GRP1/RD
      DOUBLE PRECISION X1,X2,X3,X4,A1,A2,A3,A4,A5,A6,B1,B2,B3,B4,
     1 D,DY1,DY2,DY3,DY4,AL0,AL1,AL2,AL3,DT1,DT2,DT3,DT4
      DY1=Y1
      DY2=Y2
      DY3=Y3
      DY4=Y4
      DT1=T1
      DT2=T2
      DT3=T3
      DT4=T4
      X1=DSIN(DT1)**2
      X2=DSIN(DT2)**2
      X3=DSIN(DT3)**2
      X4=DSIN(DT4)**2
      A1=X3**2*X4**2*(X4-X3)
      A2=X2**2*X4**2*(X4-X2)
      A3=X2**2*X3**2*(X3-X2)
      A4=X1**2*X4**2*(X4-X1)
      A5=X1**2*X3**2*(X3-X1)
      A6=X1**2*X2**2*(X2-X1)
      B1=X2*A1-X3*A2+X4*A3
      B2=X1*A1-X3*A4+X4*A5
      B3=X1*A2-X2*A4+X4*A6
      B4=X1*A3-X2*A5+X3*A6
      D=B1-B2+B3-B4
      AL0=(DY1*B1-DY2*B2+DY3*B3-DY4*B4)/D
      AL1=((DY2-AL0)*A1-(DY3-AL0)*A2+(DY4-AL0)*A3)/B1
      AL2=((DY3-AL0-AL1*X3)*X4**3-(DY4-AL0-AL1*X4)*X3**3)/A1
      AL3=(DY2-AL0-AL1*X2-AL2*X2**2)/X2**3
      RETURN
      END
