source: trunk/pgm97/analy.f @ 1

Last change on this file since 1 was 1, checked in by cbipsl, 18 years ago

Geisa inital import
File size: 10.9 KB
Line 
1C     CETTE OPTION DOIT VENIR APRES LA MODIFICATION PAR
2C     &GEISA PGM='TRS',UNITE='BINAIRE',MODIF='OUI' &END
3C     LE PROGRAMME ANL MODIFIE LA BANQUE DE LA FACON SUIVANTE POUR
4C     STOCKER LES FREQUENCES DES ISOTOPES DES MOLECULES :
5C     EN RECORD 1 AJOUTER APRES LL1 LA VALEUR  DE LL2=1
6C     A PARTIR DU RECORD IFIN+LL1 ECRIRE LES VALEURS SUIVANTES :
7C     NTAB,(TAB(J),J=1,NTAB),(YMOY(J),J=1,20),(YMAX(J),J=1,20),
8C     (ALF(J),J=1,20)
9C     MODE=-1  APPEL NORMAL DE ANL POUR LISTER LES FREQUENCES MOLECULES
10C              ENTRE NU1 ET NU2
11C     MODE=0   MODIFICATION DE LA BANQUE (VOIR PRECEDEMMENT) DANS CE CAS
12C              MODIF='OUI'
13C     MODE=1   LISTE DES FREQUENCES  DE LA BANQUE PAR MOLECULE SANS
14C              LECTURE DU FICHIER(OPTION PROVENANT DE PGM='INF')
15C       
16C     ANALYSE PAR MOLECULE ET PAR VARIETE ISOTOPIQUE DU CONTENU
17C     DE LA BANQUE DANS UN DOMAINE SPECTRAL DONNE
18C     HISTOGRAMME DES FREQUENCES DES RAIES POUR UNE CLASSE
19C     D'INTENSITE (K=1) OU DE NIVEAU DE BASE (K=2)
20C       
21C     NU1,NU2 LIMITES INF ET SUP DU DOMAINE SPECTRAL ETUDIE
22C     DNU PAS D'ETUDE(JUSQU'A 500 PAS SONT PREVUS)
23C     PAR DEFAUT TOUT LE FICHIER AVEC 1 PAS D'ETUDE
24C     XMIN(K),XMAX(K),PAS(K) RESPECTIVEMENT MINIMUM,MAXIMUM,LARGEUR
25C     DE LA CLASSE POUR LES VARIABLES INTENSITE OU NIVEAU DE BASE
26C     TAB(I) TABLEAU INDICE PAR LA VALEUR DES ISOTOPES ET CONTENANT LES
27C            FREQUENCES DES ISOTOPES
28C     CODE(I) TABLEAU CONTENANT LES CODES MOLECULES
29C     NN      TABLEAU CONTENANT LE NOMBRE D'ISOTOPES PAR MOLECULE
30C             SUIVI DES CODES ISOTOPES
31C     ALF(I)  MOYENNE DES 1/2 LARGEUR A MI-HAUTEUR PAR MOLECULE
32C     HIST(I,J) TABLEAU CONTENANT L'HISTOGRAMME A REPRESENTER
33C               I=NUMERO DE LA CLASSE
34C               J=1 HISTOGRAMME DES INTENSITES
35C               J=3 HAUTEUR DES CLASSES POUR LES INTENSITES
36C               J=2 HISTOGRAMME DES NIVEAUX DE BASE
37C               J=4 HAUTEUR DES CLASSES POUR LES NIVEAUX DE BASE
38C     ICAR  HAUTEUR MAX D'UNE CLASSE(45 CARACTERES)
39C     IVAL=0 ANALYSE POUR TOUTES MOLECULES ET TOUS ISOTOPES(OPTION PAR D
40C     IVAL#0 ANALYSE POUR UNE MOLECULE OU UN ISOTOPE
41C            IVAL=NUMERO DE LA MOLECULE OU DE L'ISOTOPE
42C     IRAN=1 POUR MOLECULE
43C     IRAN=2 POUR ISOTOPE
44C     MOLE=MOLECULE DEMANDEE EXEMPLE MOLE='H2O' OU 'CO2' ...
45C     ISOT=ISOTOPE DEMANDE EXEMPLE ISOT=161 OU 162 ...
46C     ANAL='OUI' ANALYSE PAR MOLECULE
47C     ANAL=AUTRE PAS D'ANALYSE
48C     HISTO='INTE' HISTOGRAMME INTENSITE (KHIST=1)
49C     HISTO='BASE' HISTOGRAMME NIVEAU DE BASE (KHIST=2)
50C     HISTO='DEUX' HISTOGRAMME INTENSITE ET NIVEAU DE BASE (KHIST=-1)
51C     HISTO=AUTRE PAS D'HISTOGRAMME (OPTION PAR DEFAUT)
52C     NBCLAS NOMBRE DE CLASSES DE L'HISTOGRAMME (PAR DEFAUT 10)
53C     NHIST NOMBRE MAXIMUM DE CLASSES(50)
54C     KANAL=0        HISTO
55C     KANAL=1        ANAL
56C     IUNI UNITE LOGIQUE CORRESPONDANT AU FICHIER SPECTRAL
57C       
58C --> SUBROUTINES APPELEES : IMPANL
59C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - *
60C  LAST MODIF:  06.05.91 PASSAGE A 75 MOLECULES DANS LES COMMON
61C       
62C - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - *
63      subroutine analy(tab,hist,xmin,xmax,ymoy,ymax,alf,qq,*)
64C       
65      logical*1 qq(1)
66      character*9 trs1,trs2
67      character*7 form,bin
68      character*4 code,ctlg,mole,blanc,base,deux
69      character*3 pgm,ianl,iext,itrs,ilst,icop,info,icre,
70     &            oui,liste,iopt,modif,trans
71      character*2 ikod,icod
72      character*1 bl
73      dimension pas(2),kode(2),alf(40)
74C  NOUVEAU TABLEAU V DE LECTURE ON PASSE DE 16 -> 29 OCTETS
75      real*8 ymoyd(40),ymaxd(40),cor,a2
76      real v(29),ymoy(40),ymax(40),xmin(40),xmax(40)
77      real nu1,nu2
78      integer hist(1),tab(1),vers
79C       
80      common/p1/ nu1,nu2,dnu,nbclas,khist,kanal,vers,mode,liste
81      common/p2/ mole(75),isot(100),nbi(75),nbt(75),iopt,ctlg,trans,
82     &           trs1,trs2
83      common/p3/ imole,iran,ival,pgm,ianl,iext,itrs,ilst,icop,info,icre
84      common/p4/ nmol,knmol,ksot,kksot,ntab,nhist,kp,lre,form,bin,modif
85      common/p5/ code(75),nn(150),nq(150),ikod(18),icod(18),jdeb(75)
86      common/inteh/ incr,pas1,pmax
87      common/entsor/iuni,juni,kuni,ient,isor,iper,nresv,maxx,blanc,oui
88C       
89      equivalence (kode(1),izot),(kode(2),imol)
90      equivalence (a,v(1)),(izot,v(15)),(imol,v(16))
91C       
92      data nbklas/500/,ibase/950/,cor/1.d50/
93C       
94C     INITIALISATION DES TABLEAUX
95C       
96c     print *,' entree dans ANL'
97      impr=0
98c     print *,' mode=',mode
99      if(mode.eq.-1) go to 6
100c     print *,' appel PGEISA'
101      call pgeisa(0.,99999.,999)
102c     print *,' retour PGEISA'
103999   read (iuni,rec=1)
104     &nu1,nu2,anu,n203,nbraie,nbmol,iecr,ifin,ll1,ll2,ll3,ll4
105      if(mode.eq.0) go to 7
106      if(mode.eq.1.and.ll2.ne.0) go to 5
107      write(isor,2000)
1082000  format(///' *inf*   cette option est uniquement valable pour la ba
109     &nque'/9x,'des donnees spectroscopiques'///)
110      go to 900
1115     continue
112      ifin=ifin+ll1
113      read (iuni,rec=ifin)
114     &ntab,(tab(j),j=1,ntab),nmol,(ymoy(j),j=1,nmol),
115     &(ymax(j),j=1,nmol),(alf(j),j=1,nmol)
116      do 12 j=1,nmol
117      ymoyd(j)=ymoy(j)*(1./cor)
118      ymaxd(j)=ymax(j)*(1./cor)
11912    continue
120      call impanl(tab,nu1,nu2,impr,ymoyd,ymaxd,alf,qq,ibase)
121      go to 900
1226     continue
123      call pgeisa(nu1,nu2,*900)
1247     continue
125      do 1 j=1,nmol
126      ymoy(j)=0.
127      ymax(j)=0.
128      ymoyd(j)=0.
129      ymaxd(j)=0.
1301     alf(j)=0.
131      nbc2=nbclas*2
132      npas=nbklas
133      do 10 j=1,npas
134      xmin(j)=1.e20
135      xmax(j)=0.
13610    continue
137      lk=ntab+nbc2
138      do 11 j=1,lk
13911    tab(j)=0
140C       
141C     CALCUL DES FREQUENCES PAR MOLECULE ET PAR ISOTOPE ET DES
142C     MOYENNES ALFMOY
143C       
144      if(dnu.ne.0.) im=(nu2-nu1)/dnu
145      if(im.gt.nbklas) nu2=nu1+float(nbklas)*dnu
146      im=0
147      anu1=nu1
148      anu2=nu1+dnu
149      if(dnu.eq.0.) anu2=nu2
150100   continue
151      call lgeisa(v,*200)
152C       
153C     TRAITEMENT PARTICULIER DES MOLECULES AYANT DES ISOTOPES DUPLIQUES
154C     MOLE= 9 10 17 18 19 25 26 28 30 31 36 42
155C       
156C     DANS LE TABLEAU TAB LES INDICES  951 A 992 SONT RESERVES POUR
157C     LES ISOTOPES DUPLIQUES (IBASE=950)
158C BB 06.05.97 cas c2h4 ( 2 isotopes dupliques le 2emendice=ibase=950
159C       
160C            H2O  CO2  O3  N2O  CO   CH4  O2   NO   SO2  NO2  NH3  PH3
161      go to (541, 541,541, 541, 541, 541, 541, 541, 509, 510, 541, 541,
162C            HNO3 OH   HF   HCL  HBR  HI   CLO  OCS  H2CO C2H6 CH3D C2H2
163     &       541, 541, 541, 541, 517, 518, 519, 541, 541, 541, 541, 541,
164C            C2H4 GEH4 HCN  C3H8 C2N2 C4H2 HC3N HOCL N2   CH3CLH2O2 H2S
165     &       525, 526, 541, 528, 541, 530, 531, 541, 541, 541, 541, 536,
166C            HCOOH COF2 SF6 C3H4  HO2  CLONO2
167     &       541,  541, 541, 541, 537,  541 ),imol
168C       
169C     SO2
170509   continue
171      if(izot.eq.626) izot=ibase+imol
172      go to 541
173C       
174C     NO2
175510   continue
176      if(izot.eq.646) izot=ibase+imol
177      go to 541
178C       
179C     HBR
180517   continue
181      if(izot.eq. 19) izot=ibase+imol
182      go to 541
183C       
184C     HI
185518   continue
186      if(izot.eq. 17) izot=ibase+imol
187      go to 541
188C       
189C     CLO
190519   continue
191      if(izot.eq. 56) izot=ibase+imol
192      go to 541
193C       
194C     C2H4
195525   continue
196      if(izot.eq.211) izot=ibase+imol
197      if(izot.eq.311) izot=ibase
198      go to 541
199C       
200C     GEH4
201526   continue
202      if(izot.eq.411) izot=ibase+imol
203      go to 541
204C       
205C     C3H8
206528   continue
207      if(izot.eq.221) izot=ibase+imol
208      go to 541
209C       
210C     C4H2
211530   continue
212      if(izot.eq.211) izot=ibase+imol
213      go to 541
214C       
215C     HC3N
216531   continue
217      if(izot.eq.124) izot=ibase+imol
218      go to 541
219C       
220C     H2S
221536   continue
222      if(izot.eq.131) izot=ibase+imol
223C     GO TO 541
224C
225C     HO2
226537   continue
227      if(izot.eq.166) izot=ibase+imol
228C     GO TO 541
229C       
230541   continue
231      if(a.le.anu2) go to 19
232      if(kanal.eq.1) go to 16
233      im=im+1
234      anu1=anu2
235      anu2=amin1(nu2,anu2+dnu)
236      go to 20
23716    continue
238      call impanl(tab,anu1,anu2,impr,ymoyd,ymaxd,alf,qq,ibase)
239      anu1=anu2
240      anu2=amin1(nu2,anu2+dnu)
241      do 17 j=1,lk
24217    tab(j)=0
243      do 18 j=1,nmol
244      ymoy(j)=0.
245      ymax(j)=0.
246      ymoyd(j)=0.
247      ymaxd(j)=0.
24818    alf(j)=0.
249      im=im+1
25019    continue
251      if(.not.qq(imol)) go to 100
252      if(kanal.eq.0) go to 20
253      tab(izot)=tab(izot)+1
254      alf(imol)=alf(imol)+v(3)
255      a2=v(2)*(1./cor)
256      ymoyd(imol)=ymoyd(imol)+a2
257Camax1 -> dmax1 ( real 8) BB 19/11/96
258      ymaxd(imol)=dmax1(ymaxd(imol),a2)
259C     ymax(imol)=amax1(ymax(imol),v(2))
260C       
261C     CALCUL DES MIN ET DES MAX DES CLASSES POUR LES INTENSITES ET LES
262C     NIVEAUX DE BASE
263C       
264      if(khist.eq.0.or.khist.eq.1) go to 100
26520    continue
266      if(ival.ne.0.and.kode(iran).ne.ival) go to 100
267      xmax(im+1)=amax1(xmax(im+1),v(4))
268      xmin(im+1)=amin1(xmin(im+1),v(4))
269      go to 100
270200   continue
271      if(kanal.eq.0) go to 24
272      call impanl(tab,anu1,anu2,impr,ymoyd,ymaxd,alf,qq,ibase)
273      if(mode.ne.0.or.modif.ne.oui) go to 900
274c     print *,' avant read rec=1   '
275      read (iuni,rec=1)
276     &nu1,nu2,anu,n203,nbraie,nbmol,iecr,ifin,ll1,ll2,ll3,ll4
277c     print *,' apres read rec=1   '
278c     print *,
279c    &nu1,nu2,anu,n203,nbraie,nbmol,iecr,ifin,ll1,ll2,ll3,ll4
280      ll2=1
281c     print *,' avant write rec=1   '
282c     print *, ifin,ll1
283      write(iuni,rec=1)
284     &nu1,nu2,anu,n203,nbraie,nbmol,iecr,ifin,ll1,ll2,ll3,ll4
285c     print *,' apres write rec=1   '
286c     print *,
287c    &nu1,nu2,anu,n203,nbraie,nbmol,iecr,ifin,ll1,ll2,ll3,ll4     
288c     print *, ifin,ll1
289      ifin=ifin+ll1
290      do 22 j=1,nmol
291      ymoy(j)=ymoyd(j)*cor
29222    ymax(j)=ymaxd(j)*cor
293c     print *, ifin
294c     write(iuni,rec=ifin)
295c    & ntab,(tab(j),j=1,ntab),nmol,(ymoy(j),j=1,nmol),
296c    &(ymax(j),j=1,nmol),(alf(j),j=1,nmol)
297c     print *,
298c    & ntab,(tab(j),j=1,ntab),nmol,(ymoy(j),j=1,nmol),
299c    &(ymax(j),j=1,nmol),(alf(j),j=1,nmol)
300      go to 900
30124    continue
302      impr=0
303      im=0
304      anu1=nu1
305      anu2=nu1+dnu
306      if(dnu.eq.0.) anu2=nu2
307      pas(   2)=(xmax(im+1)-xmin(im+1))/nbclas
308      if(pas(2).eq.0.) pas(2)=1.e+20
309      call pgeisa(nu1,nu2,*900)
310300   continue
311      call lgeisa(v,*400)
312C       
313C     CALCUL DES FREQUENCES DES CLASSES DE HIST
314C       
315      if(a.le.anu2) go to 319
316      call imph  (hist,xmin(im+1),xmax(im+1),pas,anu1,anu2,impr)
317      im=im+1
318      anu1=anu2
319      anu2=amin1(nu2,anu2+dnu)
320      pas(2)=(xmax(im+1)-xmin(im+1))/nbclas
321      if(pas(2).eq.0.) pas(2)=1.e+20
322      do 26 i=1,nbc2
323      hist(i)=0
32426    continue
325319   continue
326      if(ival.ne.0.and.kode(iran).ne.ival) go to 300
327C alog -> dlog (real 8) BB 19/11/96
328      a2=v(2)*(1./cor)
329      kt=incr+dlog10(a2)
330C     kt=incr+alog10(v(2))
331      hist(kt)=hist(kt)+1
332      if(khist.eq.1) go to 300
333      kt=(v(4)-xmin(im+1) )/pas(2) + 1.
334      kt=min0(kt,nbclas)
335      ktk=kt+nbclas
336      hist(ktk)=hist(ktk)+1
337      go to 300
338400   continue
339      if(khist.ne.0) call imph(hist,xmin(im+1),xmax(im+1),pas
340     &,anu1,anu2,impr)
341900   continue
342      return 1
343      end
Note: See TracBrowser for help on using the repository browser.