Millepede-II V04-17-03
mphistab.f90
Go to the documentation of this file.
1
2! Code converted using TO_F90 by Alan Miller
3! Date: 2012-03-16 Time: 11:07:45
4
89
90! *************************** Histograms ******************************
91
93MODULE hmpcons
94 USE mpdef
95 IMPLICIT NONE
96
97 INTEGER(mpi), PARAMETER :: numhis=16
98 INTEGER(mpi), PARAMETER :: nbin=120
99 INTEGER(mpi), PARAMETER :: nsampl=120
100
101END MODULE hmpcons
102
104MODULE hmpdata
105 USE hmpcons
106 IMPLICIT NONE
107
108 INTEGER(mpi) :: lun=7
109 INTEGER(mpi) :: inhist(nbin,numhis)
110 INTEGER(mpi) :: jnhist(5,numhis)
111 INTEGER(mpi) :: khist(numhis)=0
112 INTEGER(mpi) :: kvers(numhis)
113 REAL(mps) :: fnhist(nsampl,numhis)
114 REAL(mps) :: xl(6,numhis)
115 REAL(mpd) :: dl(2,numhis)
116 CHARACTER (LEN=60):: htext(numhis)
117
118END MODULE hmpdata
119
121SUBROUTINE hmpdef(ih,xa,xb,text)
122 USE hmpdata
123
124 IMPLICIT NONE
125
126 ! book millepede histogram, 'nbin' bins
127
128 INTEGER(mpi), INTENT(IN) :: ih
129 REAL(mps), INTENT(IN) :: xa
130 REAL(mps), INTENT(IN) :: xb
131 CHARACTER (LEN=*), INTENT(IN) :: text
132
133 SAVE
134 ! ...
135 IF(ih <= 0.OR.ih > numhis) RETURN
136 ! IF(XA.EQ.XB) RETURN
137 inhist(:,ih)=0
138 jnhist(:,ih)=0
139 xl(1,ih)=xa
140 xl(2,ih)=xb
141 xl(3,ih)=0.0
142 IF(xa /= xb) xl(3,ih)=real(nbin,mps)/(xb-xa)
143 xl(6,ih)=0.5*(xa+xb) ! center
144 IF(khist(ih) == 0) THEN
145 kvers(ih)=0
146 ELSE
147 kvers(ih)=kvers(ih)+1
148 END IF
149 khist(ih)=1 ! flt.pt. (lin)
150 htext(ih)=text
151 dl(1,ih)=0.0_mpd
152 dl(2,ih)=0.0_mpd
153 RETURN
154END SUBROUTINE hmpdef
155
157SUBROUTINE hmpldf(ih,text)
158 USE hmpdata
159
160 IMPLICIT NONE
161
162 INTEGER(mpi), INTENT(IN) :: ih
163 CHARACTER (LEN=*), INTENT(IN) :: text
164
165 SAVE
166 IF(ih <= 0.OR.ih > numhis) RETURN
167 inhist(:,ih)=0
168 jnhist(:,ih)=0
169 IF(khist(ih) == 0) THEN
170 kvers(ih)=0
171 ELSE
172 kvers(ih)=kvers(ih)+1
173 END IF
174 khist(ih)=2 ! integer log
175 htext(ih)=text
176 xl(1,ih)=0.0
177 xl(2,ih)=6.0
178 RETURN
179END SUBROUTINE hmpldf
180
182SUBROUTINE hmpent(ih,x)
183 USE hmpdata
184
185 IMPLICIT NONE
186 INTEGER(mpi) :: i
187 INTEGER(mpi) :: j
188
189 INTEGER(mpi), INTENT(IN) :: ih
190 REAL(mps), INTENT(IN) :: x
191
192 IF(ih <= 0.OR.ih > numhis) RETURN
193 IF(khist(ih) /= 1) RETURN
194 IF(jnhist(4,ih) >= 2147483647) RETURN
195 jnhist(4,ih)=jnhist(4,ih)+1 ! count
196 IF(jnhist(4,ih) <= nsampl) THEN
197 fnhist(jnhist(4,ih),ih)=x ! store value
198 IF(jnhist(4,ih) == nsampl) THEN
199 CALL hmpmak(inhist(1,ih),fnhist(1,ih),jnhist(1,ih), xl(1,ih),dl(1,ih))
200 END IF
201 RETURN
202 END IF
203 ! IF(JNHIST(1,IH)+JNHIST(2,IH)+JNHIST(3,IH).EQ.0) THEN
204 ! XL(4,IH)=X
205 ! XL(5,IH)=X
206 ! END IF
207 i=int(1.0+xl(3,ih)*(x-xl(1,ih)),mpi) ! X - Xmin
208 j=2
209 IF(i < 1) j=1
210 IF(i > nbin) j=3
211 jnhist(j,ih)=jnhist(j,ih)+1
212 xl(4,ih)=min(xl(4,ih),x)
213 xl(5,ih)=max(xl(5,ih),x)
214 IF(j /= 2) RETURN
215 inhist(i,ih)=inhist(i,ih)+1
216 dl(1,ih)=dl(1,ih)+ x-xl(6,ih)
217 dl(2,ih)=dl(2,ih)+(x-xl(6,ih))**2
218 RETURN
219END SUBROUTINE hmpent
220
222SUBROUTINE hmplnt(ih,ix)
223 USE hmpdata
224
225 IMPLICIT NONE
226 INTEGER(mpi) :: i
227 INTEGER(mpi) :: j
228
229 INTEGER(mpi), INTENT(IN) :: ih
230 INTEGER(mpi), INTENT(IN) :: ix
231
232 IF(ih <= 0.OR.ih > numhis) RETURN
233 IF(khist(ih) /= 2) RETURN
234 IF(jnhist(1,ih) >= 2147483647) RETURN
235 IF(ix <= 0) THEN
236 jnhist(1,ih)=jnhist(1,ih)+1
237 ELSE
238 IF(jnhist(4,ih) == 0) jnhist(4,ih)=ix
239 IF(jnhist(5,ih) == 0) jnhist(5,ih)=ix
240 jnhist(4,ih)=min(jnhist(4,ih),ix)
241 jnhist(5,ih)=max(jnhist(5,ih),ix)
242 i=int(1.0+20.0*log10(real(ix,mps)),mpi)
243 j=2
244 IF(i < 1) j=1
245 IF(i > nbin) j=3
246 IF(j == 2) inhist(i,ih)=inhist(i,ih)+1
247 jnhist(j,ih)=jnhist(j,ih)+1
248 END IF
249 RETURN
250END SUBROUTINE hmplnt
251
253SUBROUTINE hmprnt(ih)
254 USE hmpdata
255
256 IMPLICIT NONE
257 INTEGER(mpi) :: iha
258 INTEGER(mpi) :: ihb
259 INTEGER(mpi) :: ihc
260 INTEGER(mpi) :: j
261 INTEGER(mpi) :: nn
262 REAL(mps) :: xcent
263 REAL(mps) :: xmean
264 REAL(mps) :: xsigm
265
266 INTEGER(mpi), INTENT(IN) :: ih
267
268 IF(ih == 0) THEN
269 iha=1
270 ihb=numhis
271 ELSE
272 IF(ih <= 0.OR.ih > numhis) RETURN
273 iha=ih
274 ihb=ih
275 END IF
276 DO ihc=iha,ihb
277 IF(khist(ihc) /= 0) THEN
278 IF(khist(ihc) == 1) THEN
279 CALL hmpmak(inhist(1,ihc),fnhist(1,ihc),jnhist(1,ihc), &
280 xl(1,ihc),dl(1,ihc))
281 END IF
282 nn=jnhist(1,ihc)+jnhist(2,ihc)+jnhist(3,ihc)
283 IF(nn /= 0.OR.khist(ihc) == 3) THEN
284 WRITE(*,111)
285111 FORMAT(' ______',2('______________________________'))
286 IF(kvers(ihc) == 1) THEN
287 WRITE(*,*) 'Histogram',ihc,': ',htext(ihc)
288 ELSE
289 WRITE(*,*) 'Histogram',ihc,'/',kvers(ihc),': ',htext(ihc)
290 END IF
291 IF(khist(ihc) == 1) THEN
292 WRITE(*,*) ' Out_low inside out_high = ', (jnhist(j,ihc),j=1,3)
293 ELSE IF(khist(ihc) == 2) THEN
294 WRITE(*,*) ' 0_or_negative inside above_10^6 = ', &
295 (jnhist(j,ihc),j=1,3)
296 END IF
297 IF(khist(ihc) == 3) THEN
298 CALL pfvert(nbin,fnhist(1,ihc))
299 END IF
300 IF(jnhist(2,ihc) /= 0) THEN ! integer content
301 CALL pivert(nbin,inhist(1,ihc))
302 IF(khist(ihc) == 1) THEN
303 CALL psvert(xl(1,ihc),xl(2,ihc))
304 ELSE IF(khist(ihc) == 2) THEN
305 CALL psvert(0.0,6.0)
306 END IF
307 END IF
308 IF(khist(ihc) == 1) THEN
309 WRITE(*,*) ' Min and Max are',xl(4,ihc),xl(5,ihc)
310 IF(jnhist(2,ihc) > 1) THEN
311 xmean=real(xl(6,ihc)+dl(1,ihc)/real(jnhist(2,ihc),mps),mps)
312 xcent=0.5*(xl(1,ihc)+xl(2,ihc))
313 xsigm=real((dl(2,ihc)-dl(1,ihc)**2/real(jnhist(2,ihc),mps)),mps)
314 xsigm=sqrt(xsigm/real(jnhist(2,ihc)-1,mps))
315 WRITE(*,*) ' Mean and sigma are', xmean,' +-',xsigm
316 END IF
317 ELSE IF(khist(ihc) == 2) THEN
318 WRITE(*,*) ' Plot of log10 of entries. Min and Max are', &
319 jnhist(4,ihc),jnhist(5,ihc)
320 END IF
321 END IF
322 END IF
323 END DO
324 RETURN
325END SUBROUTINE hmprnt
326
328SUBROUTINE hmplun(lunw)
329 USE hmpdata
330
331 IMPLICIT NONE
332
333 INTEGER(mpi), INTENT(IN) :: lunw
334
335 lun=lunw
336 RETURN
337END SUBROUTINE hmplun
338
340SUBROUTINE hmpwrt(ih)
341 USE hmpdata
342
343 IMPLICIT NONE
344 INTEGER(mpi) :: iha
345 INTEGER(mpi) :: ihb
346 INTEGER(mpi) :: ihc
347 INTEGER(mpi) :: i
348 INTEGER(mpi) :: j
349 REAL(mps) :: xcent
350 REAL(mps) :: xmean
351 REAL(mps) :: xsigm
352
353 INTEGER(mpi), INTENT(IN) :: ih
354
355 IF(lun <= 0) RETURN
356 IF(ih == 0) THEN
357 iha=1
358 ihb=numhis
359 ELSE
360 IF(ih <= 0.OR.ih > numhis) RETURN
361 iha=ih
362 ihb=ih
363 END IF
364
365 DO ihc=iha,ihb ! histogram loop
366 IF(khist(ihc) /= 0) THEN
367 IF(khist(ihc) == 1) THEN
368 CALL hmpmak(inhist(1,ihc),fnhist(1,ihc),jnhist(1,ihc), &
369 xl(1,ihc),dl(1,ihc))
370 END IF
371 WRITE(lun,204) ' '
372 WRITE(lun,201) ihc,kvers(ihc),khist(ihc)
373 WRITE(lun,204) htext(ihc)
374 IF (jnhist(1,ihc)+jnhist(2,ihc)+jnhist(3,ihc) == 0 &
375 .AND.xl(1,ihc) == xl(2,ihc)) THEN
376 ! hist is empty and hist range makes no sense
377 ! - cause: hist with 'variable edges' was never filled
378 ! - workaround: make lower and upper edge of hist differ in output
379 WRITE(lun,202) nbin,xl(1,ihc)-0.001,xl(2,ihc)+0.001
380 ELSE
381 WRITE(lun,202) nbin,xl(1,ihc),xl(2,ihc)
382 END IF
383 WRITE(lun,203) (jnhist(j,ihc),j=1,3)
384 WRITE(lun,204) 'bincontent'
385 IF(khist(ihc) == 1.OR.khist(ihc) == 2) THEN
386 CALL kprint(lun,inhist(1,ihc),nbin)
387 ELSE
388 WRITE(lun,219) (fnhist(i,ihc),i=1,nbin)
389 END IF
390
391 IF(khist(ihc) == 1) THEN
392 WRITE(lun,205) xl(4,ihc),xl(5,ihc)
393 ELSE IF(khist(ihc) == 2) THEN
394 WRITE(lun,205) real(jnhist(4,ihc),mps),real(jnhist(5,ihc),mps)
395 END IF
396 IF(khist(ihc) == 1) THEN
397 IF(jnhist(2,ihc) > 1) THEN
398 xmean=real(xl(6,ihc)+dl(1,ihc)/real(jnhist(2,ihc),mps),mps)
399 xcent=0.5*(xl(1,ihc)+xl(2,ihc))
400 xsigm=real((dl(2,ihc)-dl(1,ihc)**2/real(jnhist(2,ihc),mps)),mps)
401 xsigm=sqrt(xsigm/real(jnhist(2,ihc)-1,mps))
402 WRITE(lun,206) xmean,xsigm
403 END IF
404 END IF
405 WRITE(lun,204) 'end of histogram'
406 END IF
407 END DO
408 RETURN
409
410201 FORMAT('Histogram ',i4,10x,'version ',i4,10x,'type',i2)
411202 FORMAT(10x,' bins, limits ',i4,2g15.5)
412203 FORMAT(10x,'out-low inside out-high ',3i10)
413204 FORMAT(a)
414205 FORMAT('minmax',2e15.7)
415206 FORMAT('meansigma',2e15.7)
416
417219 FORMAT(4e15.7)
418END SUBROUTINE hmpwrt
419
421SUBROUTINE hmpmak(inhist,fnhist,jnhist,xl,dl)
422 USE hmpcons
423
424 IMPLICIT NONE
425 INTEGER(mpi) :: i
426 INTEGER(mpi) :: j
427 INTEGER(mpi) :: k
428 INTEGER(mpi) :: nn
429 REAL(mps) :: x
430 REAL(mps) :: xa
431 REAL(mps) :: xb
432
433 INTEGER(mpi), INTENT(OUT) :: inhist(nbin)
434 REAL(mps), INTENT(IN) :: fnhist(nsampl)
435 INTEGER(mpi), INTENT(IN OUT) :: jnhist(5)
436 REAL(mps), INTENT(IN OUT) :: xl(6)
437 REAL(mpd), INTENT(OUT) :: dl(2)
438 REAL(mps) :: cphist(nsampl)
439
440 SAVE
441 ! ...
442 nn=jnhist(4)
443 ! WRITE(*,*) 'HMPMAK: NN,JNHIST(5)',NN,JNHIST(5)
444 IF(nn == 0.OR.jnhist(5) /= 0) RETURN
445 jnhist(5)=1
446 cphist(:)=fnhist(:)
447 CALL heapf(cphist,nn)
448 IF(xl(3) == 0.0) THEN
449 CALL bintab(cphist,nn,xa,xb)
450 xl(1)=xa
451 xl(2)=xb
452 xl(3)=0.0
453 IF(xa /= xb) xl(3)=real(nbin,mps)/(xb-xa)
454 xl(6)=0.5*(xa+xb) ! center
455 END IF
456 xl(4)=cphist( 1)
457 xl(5)=cphist(nn)
458 ! WRITE(*,*) 'XL ',XL
459 DO i=1,nn
460 inhist(i)=0
461 END DO
462 DO k=1,nn
463 x=cphist(k)
464 i=int(1.0+xl(3)*(x-xl(1)),mpi) ! X - Xmin
465 ! WRITE(*,*) 'K,I,X ',K,I,X
466 j=2
467 IF(i < 1) j=1
468 IF(i > nbin) j=3
469 jnhist(j)=jnhist(j)+1
470 IF(j == 2) THEN
471 inhist(i)=inhist(i)+1
472 dl(1)=dl(1)+ x-xl(6)
473 dl(2)=dl(2)+(x-xl(6))**2
474 END IF
475 END DO
476END SUBROUTINE hmpmak
477
479SUBROUTINE bintab(tab,n,xa,xb)
480 USE hmpcons
481
482 IMPLICIT NONE
483 REAL(mps) :: dd
484 REAL(mps) :: dx
485 INTEGER(mpi) :: i
486 INTEGER(mpi) :: iexp
487 INTEGER(mpi) :: ii
488 INTEGER(mpi) :: j
489 INTEGER(mpi) :: m1
490 INTEGER(mpi) :: m2
491 INTEGER(mpi) :: n1
492 INTEGER(mpi) :: n2
493 INTEGER(mpi) :: nch
494 REAL(mps) :: rat
495 REAL(mps) :: x1
496 REAL(mps) :: x2
497 REAL(mps) :: xx
498 ! Bin limits XA and XB from TAB(N)
499
500 INTEGER(mpi), INTENT(IN) :: n
501 REAL(mps), INTENT(IN) :: tab(n)
502 REAL(mps), INTENT(OUT) :: xa
503 REAL(mps), INTENT(OUT) :: xb
504
505 REAL(mps) :: bin(10)
506 DATA bin/1.0,1.5,2.0,3.0,4.0,5.0,8.0,10.0,15.0,20.0/
507 SAVE
508 ! ...
509
510 CALL heapf(tab,n) ! reduced statistic
511 ! WRITE(*,*) ' '
512 ! WRITE(*,*) 'Sorted ',N,(TAB(I),I=1,N)
513 IF(n < nsampl) THEN
514 x1=tab(1)
515 x2=tab(n)
516 ! WRITE(*,*) 'reduced statistic X1 X2 ',X1,X2
517 ELSE ! large statistic
518 m1=int(1.0+0.05*real(n),mpi)
519 m2=int(1.0+0.16*real(n),mpi)
520 x1=tab(m1)-4.0*(tab(m2)-tab(m1))
521 IF(x1 < 0.0.AND.tab(1) >= 0.0) x1=tab(1)
522 x2=tab(n+1-m1)+4.0*(tab(n+1-m1)-tab(n+1-m2))
523 IF(x2 > 0.0.AND.tab(n) <= 0.0) x2=tab(n)
524 ! WRITE(*,*) 'large statistic ',X1,X2
525 ! WRITE(*,*) 'min und max ',TAB(1),TAB(N)
526 IF(x1*tab(1) <= 0.0) x1=0.0
527 IF(x2*tab(n) <= 0.0) x2=0.0
528 ! WRITE(*,*) 'large statistic zero ',X1,X2
529 IF(x1*x2 < 0.0.AND.min(-x1,x2) > 0.6*max(-x1,x2)) THEN
530 xx=max(-x1,x2) ! symmetry
531 x1=-xx
532 x2=+xx
533 ELSE IF(x1*x2 > 0.0.AND. & ! include zero ?
534 abs(min(x1,x2)) < 0.4*abs(max(x1,x2))) THEN
535 IF(x1 < 0.0) THEN
536 x2=0.0
537 ELSE
538 x1=0.0
539 END IF
540 END IF
541 ! WRITE(*,*) 'large statistic ',X1,X2
542 END IF
543 IF(x1 == x2) THEN
544 x1=x1-1.0
545 x2=x2+1.0
546 END IF
547 dx=x2-x1
548 ! WRITE(*,*) 'X1,X2,DX ',N,X1,X2,DX
549 rat=0.0
550 ii=1
551 DO j=1,11
552 i=j
553 IF(j == 11) i=ii
554 iexp=int(101.0+log10(dx)-log10(6.0*bin(i)),mpi)
555 iexp=iexp-100
556 dd=bin(i)*10.0**iexp
557
558 n1=int(abs(x1)/dd,mpi)
559 IF(x1 < 0.0) n1=-n1
560 IF(real(n1,mps)*dd > x1) n1=n1-1
561 ! WRITE(*,*) 'Bin ',I,N1,N1*DD,X1
562
563 n2=int(abs(x2)/dd,mpi)
564 IF(x2 < 0.0) n2=-n2
565 IF(real(n2,mps)*dd < x2) n2=n2+1
566 ! WRITE(*,*) 'Bin ',I,N1,N2,N2*DD,X2
56710 IF(n2-n1 < 6) THEN
568 nch=0 ! number of changes
569 IF(n1 /= 0) THEN
570 n1=n1-1
571 nch=nch+1
572 END IF
573 IF(n2-n1 < 6.AND.n2 /= 0) THEN
574 n2=n2+1
575 nch=nch+1
576 ENDIF
577 IF (nch > 0) GO TO 10
578 ! avoid infinite loop
579 print *, ' BINTAB: break infinite loop ', n1, n2, n, x1, x2, dd
580 n1=-3
581 n2=3
582 END IF
583 ! WRITE(*,*) 'corrected N1 N2 ',N1,N2
584 xa=sign(real(n1,mps)*dd,x1)
585 xb=sign(real(n2,mps)*dd,x2)
586 ! WRITE(*,*) J,' resulting limits XA XB ',XA,XB
587 IF((x2-x1)/(xb-xa) > rat) THEN
588 ii=i
589 rat=(x2-x1)/(xb-xa)
590 END IF
591 END DO
592! WRITE(*,*) J,' resulting limits XA XB ',XA,XB
593END SUBROUTINE bintab
594
596SUBROUTINE kprint(lun,list,n)
597 USE mpdef
598
599 IMPLICIT NONE
600 INTEGER(mpi) :: i
601 INTEGER(mpi) :: ia
602 INTEGER(mpi) :: ib
603 INTEGER(mpi) :: k
604 INTEGER(mpi) :: ln
605 INTEGER(mpi) :: lp
606 INTEGER(mpi) :: np
607 ! print integer array LIST(N)
608
609 INTEGER(mpi), INTENT(IN OUT) :: lun
610 INTEGER(mpi), INTENT(IN) :: n
611 INTEGER(mpi), INTENT(IN) :: list(n)
612 INTEGER(mpi) :: li(7)
613 DATA li/2,3,4,6,8,9,12/ ! number of characters
614 SAVE
615 ! ...
616 ib=0
61710 ia=ib+1
618 IF(ia > n) RETURN
619 DO k=1,7
620 np=72/li(k)
621 ib=min(ia-1+np,n)
622 IF(k <= 6) THEN
623 lp=10**(li(k)-1)-1 ! maximum positive
624 ln=-lp/10 ! minimum negative
625 DO i=ia,ib
626 IF(list(i) > lp.OR.list(i) < ln) GO TO 20
627 END DO
628 END IF
629 IF(k == 1) THEN
630 WRITE(lun,101) (list(i),i=ia,ib)
631 ELSE IF(k == 2) THEN
632 WRITE(lun,102) (list(i),i=ia,ib)
633 ELSE IF(k == 3) THEN
634 WRITE(lun,103) (list(i),i=ia,ib)
635 ELSE IF(k == 4) THEN
636 WRITE(lun,104) (list(i),i=ia,ib)
637 ELSE IF(k == 5) THEN
638 WRITE(lun,105) (list(i),i=ia,ib)
639 ELSE IF(k == 6) THEN
640 WRITE(lun,106) (list(i),i=ia,ib)
641 ELSE IF(k == 7) THEN
642 WRITE(lun,107) (list(i),i=ia,ib)
643 END IF
644 GO TO 10
64520 CONTINUE
646 END DO
647101 FORMAT(36i2)
648102 FORMAT(24i3)
649103 FORMAT(18i4)
650104 FORMAT(12i6)
651105 FORMAT( 9i8)
652106 FORMAT( 8i9)
653107 FORMAT( 6i12)
654END SUBROUTINE kprint
655
656! ***************************** XY data ****************************
657
659MODULE gmpcons
660 USE mpdef
661 IMPLICIT NONE
662
663 INTEGER(mpi), PARAMETER :: numgxy=10
664
665END MODULE gmpcons
666
668MODULE gmpdata
669 USE gmpcons
670 IMPLICIT NONE
671
672 INTEGER(mpi) :: lun=7
673 INTEGER(mpi) :: nlimit=500
674 INTEGER(mpi) :: nstr(numgxy)=0
675 INTEGER(mpi) :: igtp(numgxy)
676 ! ITYP = 1 X,Y as dots
677 ! = 2 X,Y as line
678 ! = 3 X,Y as line and dots
679 ! = 4 X,Y, DX,DY symbols
680 ! = 5 mean, sigma of Y vs X (GMPMS)
681 INTEGER(mpi) :: lvers(numgxy)
682 INTEGER(mpi) :: nst(3,numgxy)
683
684 INTEGER(mpi) :: jflc(5,numgxy)
685 INTEGER(mpi) :: kflc(5,numgxy)
686 ! JFLC(1,.) = first used index
687 ! JFLC(2,.) = last used index
688 ! JFLC(3,.) = counter of used places
689 ! JFLC(4,.) = counter of ignored
690 ! JFLC(5,.) = limit for JFLC(3)
691 REAL(mps) :: xyplws(10,numgxy)
692 REAL(mps) :: y1(numgxy)
693 REAL(mps), DIMENSION(:,:), ALLOCATABLE :: array
694 REAL(mps), DIMENSION(:,:), ALLOCATABLE :: array4
695 REAL(mps), DIMENSION(:), ALLOCATABLE :: array1
696 CHARACTER (LEN=60):: gtext(numgxy)
697
698END MODULE gmpdata
699
701SUBROUTINE gmpdef(ig,ityp,text)
702 USE gmpdata
703
704 IMPLICIT NONE
705 INTEGER(mpi) :: i
706 INTEGER(mpi) :: j
707
708 INTEGER(mpi), INTENT(IN) :: ig
709 INTEGER(mpi), INTENT(IN) :: ityp
710 CHARACTER (LEN=*), INTENT(IN) :: text
711
712 LOGICAL:: start
713 SAVE
714 DATA start/.true./
715 ! ...
716 IF(start) THEN
717 start=.false.
718 ! dummy call to increase nlimit ?
719 if(ig == 0) nlimit = max(nlimit, ityp)
720 CALL stmars(nlimit*numgxy) ! initialize storage
721 DO i=1,numgxy
722 DO j=1,5
723 jflc(j,i)=0
724 kflc(j,i)=0
725 END DO
726 END DO
727 ALLOCATE (array(2,nlimit))
728 ALLOCATE (array4(4,(nlimit+1)/2))
729 ALLOCATE (array1(nlimit*2))
730 END IF
731
732 IF(ig < 1.OR.ig > numgxy) RETURN
733 IF(ityp < 1.OR.ityp > 5) RETURN
734 IF(nstr(ig) == 0) THEN
735 lvers(ig)=0
736 ELSE
737 lvers(ig)=lvers(ig)+1
738 END IF
739 nstr(ig)=1 ! by GF
740 ! remove stored elements
741 IF(jflc(1,ig) /= 0) CALL stmarm(jflc(1,ig))
742 IF(kflc(1,ig) /= 0) CALL stmarm(kflc(1,ig))
743 igtp(ig)=ityp
744 gtext(ig)=text
745 DO j=1,5
746 jflc(j,ig)=0
747 END DO
748 jflc(5,ig)=nlimit
749 IF(ityp == 5) THEN
750 DO j=1,5
751 kflc(j,ig)=0
752 END DO
753 jflc(5,ig)=128 ! maximum of 128 values
754 kflc(5,ig)=nlimit
755 nst(1,ig)=0
756 nst(2,ig)=0
757 nst(3,ig)=1
758 DO j=1,10
759 xyplws(j,ig)=0.0
760 END DO
761 END IF
762 RETURN
763END SUBROUTINE gmpdef
764
766SUBROUTINE gmpxy(ig,x,y)
767 USE gmpdata
768
769 IMPLICIT NONE
770 INTEGER(mpi), INTENT(IN) :: ig
771 REAL(mps), INTENT(IN) :: x
772 REAL(mps), INTENT(IN) :: y
773
774 IF(ig < 1.OR.ig > numgxy) RETURN ! check argument IG
775 IF(igtp(ig) < 1.OR.igtp(ig) > 3) RETURN ! check type
776 CALL stmapr(jflc(1,ig),x,y)
777 RETURN
778END SUBROUTINE gmpxy
779
781SUBROUTINE gmpxyd(ig,x,y,dx,dy)
782 USE gmpdata
783
784 IMPLICIT NONE
785 REAL(mps) :: four(4)
786
787 INTEGER(mpi), INTENT(IN) :: ig
788 REAL(mps), INTENT(IN) :: x
789 REAL(mps), INTENT(IN) :: y
790 REAL(mps), INTENT(IN) :: dx
791 REAL(mps), INTENT(IN) :: dy
792
793 IF(ig < 1.OR.ig > numgxy) RETURN ! check argument IG
794 IF(igtp(ig) /= 4) RETURN
795 four(1)=x
796 four(2)=y
797 four(3)=dx
798 four(4)=dy
799 CALL stmadp(jflc(1,ig),four)
800 RETURN
801END SUBROUTINE gmpxyd
802
804SUBROUTINE gmpms(ig,x,y)
805 USE gmpdata
806
807 IMPLICIT NONE
808 INTEGER(mpi) :: i
809 INTEGER(mpi) :: n
810
811 INTEGER(mpi), INTENT(IN) :: ig
812 REAL(mps), INTENT(IN) :: x
813 REAL(mps), INTENT(IN) :: y
814
815 ! mean sigma from Y, as a function of X
816 ! WRITE(*,*) 'GMPMS ',IG,X,Y
817
818 IF(ig < 1.OR.ig > numgxy) RETURN ! check argument IG
819 IF(igtp(ig) /= 5) RETURN
820
821 xyplws(10,ig)=x ! last X coordinate
822 IF(nst(1,ig) == 0) THEN
823 y1(ig)=y
824 nst(1,ig)=1
825 IF(kflc(3,ig) == 0) xyplws(9,ig)=x ! start coordinate
826 ELSE
827 nst(1,ig)=0
828 CALL stmapr(kflc(1,ig),y1(ig),y) ! store pair
829 IF(kflc(3,ig) >= kflc(5,ig)) THEN
830 CALL stmacp1(kflc(1,ig),array1,n) ! get data
831 IF(kflc(1,ig) /= 0) CALL stmarm(kflc(1,ig)) ! remove data
832 CALL rmesig(array1,n,xyplws(2,ig),xyplws(4,ig))
833 nst(2,ig)=nst(2,ig)+1
834 IF(nst(2,ig) == 1) xyplws(7,ig)=xyplws(9,ig)
835 xyplws(8,ig)=x ! end coordinate
836 xyplws(5,ig)=xyplws(5,ig)+xyplws(2,ig)
837 xyplws(6,ig)=xyplws(6,ig)+xyplws(4,ig)
838 IF(nst(2,ig) == nst(3,ig)) THEN
839 xyplws(1,ig)=0.5*(xyplws(7,ig)+xyplws(8,ig))
840 xyplws(2,ig)=xyplws(5,ig)/real(nst(3,ig),mps)
841 xyplws(3,ig)=0.5*(xyplws(8,ig)-xyplws(7,ig))
842 xyplws(4,ig)=xyplws(6,ig)/real(nst(3,ig),mps)
843 xyplws(5,ig)=0.0
844 xyplws(6,ig)=0.0
845 nst(2,ig)=0
846 CALL stmadp(jflc(1,ig),xyplws(1,ig))
847 IF(jflc(3,ig) >= jflc(5,ig)) THEN
848 CALL stmacp4(jflc(1,ig),array4,n) ! get data
849 IF(jflc(1,ig) /= 0) CALL stmarm(jflc(1,ig)) ! remove data
850 DO i=1,n,2 ! average
851 xyplws(7,ig)=array4(1,i )-array4(3, i)
852 xyplws(8,ig)=array4(1,i+1)+array4(3,i+1)
853 xyplws(1,ig)=0.5*(xyplws(7,ig)+xyplws(8,ig))
854 xyplws(2,ig)=0.5*(array4(2,i)+array4(2,i+1))
855 xyplws(3,ig)=0.5*(xyplws(8,ig)-xyplws(7,ig))
856 xyplws(4,ig)=0.5*(array4(4,i)+array4(4,i+1))
857 CALL stmadp(jflc(1,ig),xyplws(1,ig))
858 END DO
859 nst(3,ig)=nst(3,ig)+nst(3,ig)
860 END IF
861 END IF
862 END IF
863 END IF
864 RETURN
865END SUBROUTINE gmpms
866
868SUBROUTINE gmprnt(ig)
869 USE gmpdata
870
871 IMPLICIT NONE
872 INTEGER(mpi) :: iga
873 INTEGER(mpi) :: igb
874 INTEGER(mpi) :: igc
875 INTEGER(mpi) :: j
876 INTEGER(mpi) :: n
877 INTEGER(mpi) :: na
878 REAL(mps) :: wght
879
880 INTEGER(mpi), INTENT(IN) :: ig
881
882 IF(ig == 0) THEN
883 iga=1
884 igb=numgxy
885 ELSE
886 IF(ig <= 0.OR.ig > numgxy) RETURN
887 iga=ig
888 igb=ig
889 END IF
890 DO igc=iga,igb
891
892 IF(igtp(igc) >= 1.AND.igtp(igc) <= 3) THEN
893 WRITE(*,*) ' '
894 WRITE(*,*) 'Store ',igc,': ',gtext(igc)
895 IF(jflc(4,igc) == 0) THEN
896 WRITE(*,*) ' stored n-tuples: ',jflc(3,igc)
897 ELSE
898 WRITE(*,*) ' stored n-tuples, not-stored n-tuples: ', &
899 jflc(3,igc),', ',jflc(4,igc)
900 END IF
901
902 CALL stmacp(jflc(1,igc),array,na) ! get all data
903
904 DO n=1,na
905 WRITE(*,102) n, array(1,n),array(2,n)
906 END DO
907
908 ELSE IF(igtp(igc) == 4) THEN
909
910 WRITE(*,*) ' '
911 WRITE(*,*) 'Store ',igc,': ',gtext(igc)
912 IF(jflc(4,igc) == 0) THEN
913 WRITE(*,*) ' stored n-tuples: ',jflc(3,igc)
914 ELSE
915 WRITE(*,*) ' stored n-tuples, not-stored n-tuples: ', &
916 jflc(3,igc),', ',jflc(4,igc)
917 END IF
918
919 CALL stmacp4(jflc(1,igc),array4,na) ! get all data
920
921 DO n=1,na
922 WRITE(*,102) n,(array4(j,n),j=1,4)
923 END DO
924
925 ELSE IF(igtp(igc) == 5) THEN
926
927 CALL stmacp1(kflc(1,igc),array1,n) ! get data
928 IF(kflc(1,igc) /= 0) CALL stmarm(kflc(1,igc)) ! remove data
929 IF(nst(1,igc) == 1) THEN
930 n=n+1
931 array1(n)=y1(igc)
932 nst(1,igc)=0 ! reset
933 END IF
934 IF(n /= 0) THEN
935 xyplws(7,igc)=xyplws( 9,igc)
936 xyplws(8,igc)=xyplws(10,igc)
937 CALL rmesig(array1,n,xyplws(2,igc),xyplws(4,igc))
938 wght=real(n,mps)/real(nst(3,igc)*kflc(5,igc),mps)
939 xyplws(5,igc)=xyplws(5,igc)+xyplws(2,igc)*wght
940 xyplws(6,igc)=xyplws(6,igc)+xyplws(4,igc)*wght
941 xyplws(2,igc)=xyplws(5,igc)/(real(nst(2,igc),mps)+wght)
942 xyplws(4,igc)=xyplws(6,igc)/(real(nst(2,igc),mps)+wght)
943 xyplws(1,igc)=0.5*(xyplws(7,igc)+xyplws(8,igc))
944 xyplws(3,igc)=0.5*(xyplws(8,igc)-xyplws(7,igc))
945 CALL stmadp(jflc(1,igc),xyplws(1,igc))
946 END IF
947
948 WRITE(*,*) ' '
949 WRITE(*,*) 'Store ',igc,': ',gtext(igc)
950 IF(jflc(4,igc) == 0) THEN
951 WRITE(*,*) ' stored n-tuples: ',jflc(3,igc)
952 ELSE
953 WRITE(*,*) ' stored n-tuples, not-stored n-tuples: ', &
954 jflc(3,igc),', ',jflc(4,igc)
955 END IF
956
957 CALL stmacp4(jflc(1,igc),array4,na) ! get all data
958 DO n=1,na
959 WRITE(*,102) n,(array4(j,n),j=1,4)
960 END DO
961 END IF
962 END DO
963 RETURN
964102 FORMAT(i12,4g15.7)
965 ! 103 FORMAT(' Index ___X___ ___Y___ '/
966 ! + ' ----- -------------- --------------')
967 ! 104 FORMAT(' Index ___X___ ___Y___ ',
968 ! + ' ___DX__ ___DY__ '/
969 ! + ' ----- -------------- --------------',
970 ! + ' -------------- --------------')
971END SUBROUTINE gmprnt
972
974SUBROUTINE gmplun(lunw)
975 USE gmpdata
976
977 IMPLICIT NONE
978
979 INTEGER(mpi), INTENT(IN) :: lunw
980
981 lun=lunw
982 RETURN
983END SUBROUTINE gmplun
984
986SUBROUTINE gmpwrt(ig)
987 USE gmpdata
988
989 IMPLICIT NONE
990 INTEGER(mpi) :: iga
991 INTEGER(mpi) :: igb
992 INTEGER(mpi) :: igc
993 INTEGER(mpi) :: j
994 INTEGER(mpi) :: n
995 INTEGER(mpi) :: na
996 REAL(mps) :: wght
997
998 INTEGER(mpi), INTENT(IN) :: ig
999
1000 IF(lun <= 0) RETURN
1001 IF(ig == 0) THEN
1002 iga=1
1003 igb=numgxy
1004 ELSE
1005 IF(ig <= 0.OR.ig > numgxy) RETURN
1006 iga=ig
1007 igb=ig
1008 END IF
1009 DO igc=iga,igb
1010 IF(igtp(igc) == 5) THEN
1011
1012 CALL stmacp1(kflc(1,igc),array1,n) ! get data
1013 IF(kflc(1,igc) /= 0) CALL stmarm(kflc(1,igc)) ! remove data
1014 IF(nst(1,igc) == 1) THEN
1015 n=n+1
1016 array1(n)=y1(igc)
1017 nst(1,igc)=0 ! reset
1018 END IF
1019 IF(n /= 0) THEN
1020 xyplws(7,igc)=xyplws( 9,igc)
1021 xyplws(8,igc)=xyplws(10,igc)
1022 CALL rmesig(array1,n,xyplws(2,igc),xyplws(4,igc))
1023 wght=real(n,mps)/real(nst(3,igc)*kflc(5,igc),mps)
1024 xyplws(5,igc)=xyplws(5,igc)+xyplws(2,igc)*wght
1025 xyplws(6,igc)=xyplws(6,igc)+xyplws(4,igc)*wght
1026 xyplws(2,igc)=xyplws(5,igc)/(real(nst(2,igc),mps)+wght)
1027 xyplws(4,igc)=xyplws(6,igc)/(real(nst(2,igc),mps)+wght)
1028 xyplws(1,igc)=0.5*(xyplws(7,igc)+xyplws(8,igc))
1029 xyplws(3,igc)=0.5*(xyplws(8,igc)-xyplws(7,igc))
1030 CALL stmadp(jflc(1,igc),xyplws(1,igc))
1031 END IF
1032
1033 END IF
1034 IF(jflc(3,igc)+jflc(4,igc) /= 0) THEN
1035 WRITE(lun,204) ' '
1036 WRITE(lun,201) igc,lvers(igc),igtp(igc)
1037 WRITE(lun,204) gtext(igc)
1038 WRITE(lun,203) jflc(3,igc)+jflc(4,igc)
1039 IF(igtp(igc) >= 1.AND.igtp(igc) <= 3) THEN
1040 CALL stmacp(jflc(1,igc),array,na) ! get all data
1041 WRITE(lun,204) 'x-y'
1042 DO n=1,na
1043 WRITE(lun,205) array(1,n),array(2,n)
1044 END DO
1045 ELSE IF(igtp(igc) == 4.OR.igtp(igc) == 5) THEN
1046 CALL stmacp4(jflc(1,igc),array4,na) ! get all data
1047 WRITE(lun,204) 'x-y-dx-dy'
1048 DO n=1,na
1049 WRITE(lun,205) (array4(j,n),j=1,4)
1050 END DO
1051 END IF
1052 WRITE(lun,204) 'end of xy-data'
1053 END IF
1054 END DO
1055 RETURN
1056201 FORMAT('XY-Data ',i4,10x,'version ',i4,10x,'type',i2)
1057203 FORMAT(10x,'stored not-stored ',2i10)
1058204 FORMAT(a)
1059205 FORMAT(3x,4g15.7)
1060END SUBROUTINE gmpwrt
1061
1063MODULE stmamod
1064 USE mpdef
1065 IMPLICIT NONE
1066
1067 REAL(mps), DIMENSION(:,:), ALLOCATABLE :: tk ! pair storage for data pairs
1068 INTEGER(mpi), DIMENSION(:), ALLOCATABLE :: next ! pointer
1069 INTEGER(mpi) :: iflc1 ! first and last index of free pairs
1070 INTEGER(mpi) :: iflc2 ! first and last index of free pairs
1071
1072END MODULE stmamod
1073
1075SUBROUTINE stmars(ndim)
1076 USE stmamod
1077
1078 IMPLICIT NONE
1079 INTEGER(mpi) :: i
1080
1081 INTEGER(mpi), INTENT(IN) :: ndim
1082
1083 SAVE
1084
1085 !INTEGER(mpi) :: jflc(5) ! user array
1086 ! JFLC(1) = first used index
1087 ! JFLC(2) = last used index
1088 ! JFLC(3) = counter of used places
1089 ! JFLC(4) = counter of ignored
1090 ! JFLC(5) = limit for JFLC(3)
1091 ! ...
1092 !print *, ' stmars ndim ', ndim
1093 ALLOCATE (tk(2,ndim))
1094 ALLOCATE (next(ndim))
1095
1096 DO i=1,ndim
1097 next(i)=i+1 ! pointer to next free location
1098 tk(1,i)=0.0 ! reset
1099 tk(2,i)=0.0
1100 END DO
1101 next(ndim)=0 ! ... and end pointer
1102 iflc1=1 ! index first free pair
1103 iflc2=ndim ! index last free pair
1104 RETURN
1105END SUBROUTINE stmars ! init/
1106
1108SUBROUTINE stmapr(jflc,x,y)
1109 USE stmamod
1110
1111 IMPLICIT NONE
1112 INTEGER(mpi) :: ifre
1113
1114 INTEGER(mpi), INTENT(INOUT) :: jflc(5)
1115 REAL(mps), INTENT(IN) :: x
1116 REAL(mps), INTENT(IN) :: y
1117
1118 ifre=iflc1 ! index of free place
1119 IF(ifre == 0.OR.jflc(3) >= jflc(5)) THEN ! overflow
1120 jflc(4)=jflc(4)+1
1121 ELSE
1122 iflc1=next(ifre) ! pointer to new free location
1123 IF(jflc(1) == 0) THEN ! first item
1124 jflc(1)=ifre
1125 ELSE
1126 next(jflc(2))=ifre
1127 END IF
1128 next(ifre)=0
1129 jflc(2)=ifre ! last index
1130 jflc(3)=jflc(3)+1 ! counter
1131 tk(1,ifre)=x
1132 tk(2,ifre)=y
1133 END IF
1134 RETURN
1135END SUBROUTINE stmapr
1136
1138SUBROUTINE stmadp(jflc,four)
1139 USE stmamod
1140
1141 IMPLICIT NONE
1142 INTEGER(mpi) :: ifrea
1143 INTEGER(mpi) :: ifreb
1144
1145 INTEGER(mpi), INTENT(INOUT) :: jflc(5)
1146 REAL(mps), INTENT(IN) :: four(4)
1147
1148 ifrea=iflc1 ! index of 1. free place
1149 IF(ifrea == 0) THEN ! overflow
1150 jflc(4)=jflc(4)+1
1151 ELSE
1152 ifreb=next(iflc1) ! index of 2. free place
1153 IF(ifreb == 0.OR.jflc(3) >= 2*jflc(5)) THEN ! overflow
1154 jflc(4)=jflc(4)+1
1155 ELSE
1156 iflc1=next(ifreb) ! pointer to new free location
1157 IF(jflc(1) == 0) THEN ! first item
1158 jflc(1)=ifrea
1159 ELSE
1160 next(jflc(2))=ifrea
1161 END IF
1162 next(ifreb)=0
1163 jflc(2)=ifreb ! last index
1164 jflc(3)=jflc(3)+1 ! counter
1165 tk(1,ifrea)=four(1)
1166 tk(2,ifrea)=four(2)
1167 tk(1,ifreb)=four(3)
1168 tk(2,ifreb)=four(4)
1169 END IF
1170 END IF
1171 RETURN
1172END SUBROUTINE stmadp
1173
1175SUBROUTINE stmacp(jflc,array,n)
1176 USE stmamod
1177
1178 IMPLICIT NONE
1179 INTEGER(mpi) :: ind
1180
1181 INTEGER(mpi), INTENT(IN) :: jflc(5)
1182 INTEGER(mpi), INTENT(OUT) :: n
1183 REAL(mps), INTENT(OUT) :: array(2,*)
1184
1185 n=0
1186 ind=jflc(1)
118710 IF(ind == 0) RETURN
1188 n=n+1
1189 array(1,n)=tk(1,ind)
1190 array(2,n)=tk(2,ind)
1191 ind=next(ind)
1192 GO TO 10
1193END SUBROUTINE stmacp
1194
1196SUBROUTINE stmacp4(jflc,array,n)
1197 USE stmamod
1198
1199 IMPLICIT NONE
1200 INTEGER(mpi) :: ind1
1201 INTEGER(mpi) :: ind2
1202
1203 INTEGER(mpi), INTENT(IN) :: jflc(5)
1204 INTEGER(mpi), INTENT(OUT) :: n
1205 REAL(mps), INTENT(OUT) :: array(4,*)
1206
1207 n=0
1208 ind1=jflc(1)
120910 IF(ind1 == 0) RETURN
1210 ind2=next(ind1) ! 2nd pair
1211 IF(ind2 == 0) RETURN
1212 n=n+1
1213 array(1,n)=tk(1,ind1)
1214 array(2,n)=tk(2,ind1)
1215 array(3,n)=tk(1,ind2)
1216 array(4,n)=tk(2,ind2)
1217 ind1=next(ind2)
1218 GO TO 10
1219END SUBROUTINE stmacp4
1220
1222SUBROUTINE stmacp1(jflc,array,n)
1223 USE stmamod
1224
1225 IMPLICIT NONE
1226 INTEGER(mpi) :: ind
1227
1228 INTEGER(mpi), INTENT(IN) :: jflc(5)
1229 INTEGER(mpi), INTENT(OUT) :: n
1230 REAL(mps), INTENT(OUT) :: array(*)
1231
1232 n=0
1233 ind=jflc(1)
123410 IF(ind == 0) RETURN
1235 n=n+1
1236 array(n)=tk(1,ind)
1237 n=n+1
1238 array(n)=tk(2,ind)
1239 ind=next(ind)
1240 GO TO 10
1241END SUBROUTINE stmacp1
1242
1244SUBROUTINE stmarm(jflc)
1245 USE stmamod
1246
1247 IMPLICIT NONE
1248 INTEGER(mpi) :: j
1249
1250 INTEGER(mpi), INTENT(INOUT) :: jflc(5)
1251
1252 next(iflc2)=jflc(1) ! connect to free space
1253 iflc2=jflc(2) ! new last free index
1254 DO j=1,4
1255 jflc(j)=0
1256 END DO
1257END SUBROUTINE stmarm ! init/
1258
1260SUBROUTINE rmesig(x,n,xloc,xsca)
1261 USE mpdef
1262
1263 IMPLICIT NONE
1264 INTEGER(mpi) :: i
1265 ! robust determination of location and scale parameter,
1266 ! for Gaussian data: location=mean and scale=standard deviation
1267 ! XLOC = median of X_i (N values in array X(N))
1268 ! XCSA = median of | X_i - XLOC |, times 1.4826
1269
1270 INTEGER(mpi), INTENT(IN) :: n
1271 REAL(mps), INTENT(IN OUT) :: x(n) ! input array, modified
1272 REAL(mps), INTENT(OUT) :: xloc
1273 REAL(mps), INTENT(OUT) :: xsca
1274 SAVE
1275 ! ...
1276 xloc=0.0
1277 xsca=0.0
1278 IF(n <= 0) RETURN
1279 CALL heapf(x,n) ! sort
1280 xloc=0.5*(x((n+1)/2)+x((n+2)/2)) ! location
1281 DO i=1,n
1282 x(i)=abs(x(i)-xloc)
1283 END DO
1284 CALL heapf(x,n) ! sort
1285 xsca=1.4826*0.5*(x((n+1)/2)+x((n+2)/2)) ! dispersion
1286END SUBROUTINE rmesig
1287
1288
1289
hmplun
subroutine hmplun(lunw)
unit for output
Definition: mphistab.f90:329
hmpmak
subroutine hmpmak(inhist, fnhist, jnhist, xl, dl)
hist scale from data
Definition: mphistab.f90:422
gmpdef
subroutine gmpdef(ig, ityp, text)
book, reset XY storage
Definition: mphistab.f90:702
gmpxy
subroutine gmpxy(ig, x, y)
add (X,Y) pair
Definition: mphistab.f90:767
bintab
subroutine bintab(tab, n, xa, xb)
hist scale from data
Definition: mphistab.f90:480
hmpdef
subroutine hmpdef(ih, xa, xb, text)
book, reset histogram
Definition: mphistab.f90:122
rmesig
subroutine rmesig(x, n, xloc, xsca)
robust mean and sigma
Definition: mphistab.f90:1261
stmars
subroutine stmars(ndim)
init/reset storage
Definition: mphistab.f90:1076
gmplun
subroutine gmplun(lunw)
unit for output
Definition: mphistab.f90:975
gmpxyd
subroutine gmpxyd(ig, x, y, dx, dy)
add (X,Y,DX,DY)
Definition: mphistab.f90:782
stmadp
subroutine stmadp(jflc, four)
store double pair
Definition: mphistab.f90:1139
kprint
subroutine kprint(lun, list, n)
print integer array
Definition: mphistab.f90:597
stmapr
subroutine stmapr(jflc, x, y)
store pair (X,Y)
Definition: mphistab.f90:1109
hmpwrt
subroutine hmpwrt(ih)
write histogram text file
Definition: mphistab.f90:341
gmpwrt
subroutine gmpwrt(ig)
write XY text file
Definition: mphistab.f90:987
hmpldf
subroutine hmpldf(ih, text)
book, reset log histogram
Definition: mphistab.f90:158
stmacp
subroutine stmacp(jflc, array, n)
copy (cp) all pairs to array
Definition: mphistab.f90:1176
gmprnt
subroutine gmprnt(ig)
print XY data
Definition: mphistab.f90:869
hmpent
subroutine hmpent(ih, x)
entry flt.pt.
Definition: mphistab.f90:183
hmplnt
subroutine hmplnt(ih, ix)
entry integer
Definition: mphistab.f90:223
stmarm
subroutine stmarm(jflc)
remove (rm) stored pairs
Definition: mphistab.f90:1245
gmpms
subroutine gmpms(ig, x, y)
mean sigma(X) from Y
Definition: mphistab.f90:805
stmacp1
subroutine stmacp1(jflc, array, n)
copy (cp) all pairs to array1
Definition: mphistab.f90:1223
hmprnt
subroutine hmprnt(ih)
print, content vert
Definition: mphistab.f90:254
stmacp4
subroutine stmacp4(jflc, array, n)
copy (cp) all pairs to array4
Definition: mphistab.f90:1197
heapf
subroutine heapf(a, n)
Heap sort direct (real).
Definition: mpnum.f90:1660
gmpcons
XY constants.
Definition: mphistab.f90:659
gmpcons::numgxy
integer(mpi), parameter numgxy
number of XY data plots
Definition: mphistab.f90:663
gmpdata
Histogram data.
Definition: mphistab.f90:668
gmpdata::kflc
integer(mpi), dimension(5, numgxy) kflc
meta data
Definition: mphistab.f90:685
gmpdata::lun
integer(mpi) lun
unit for output
Definition: mphistab.f90:672
gmpdata::gtext
character(len=60), dimension(numgxy) gtext
text
Definition: mphistab.f90:696
gmpdata::array
real(mps), dimension(:,:), allocatable array
X,Y.
Definition: mphistab.f90:693
gmpdata::array4
real(mps), dimension(:,:), allocatable array4
X,Y,DX,DY.
Definition: mphistab.f90:694
gmpdata::array1
real(mps), dimension(:), allocatable array1
Y(X)
Definition: mphistab.f90:695
gmpdata::nlimit
integer(mpi) nlimit
max.
Definition: mphistab.f90:673
gmpdata::jflc
integer(mpi), dimension(5, numgxy) jflc
meta data
Definition: mphistab.f90:684
gmpdata::y1
real(mps), dimension(numgxy) y1
first Y (as X) for GMPMS
Definition: mphistab.f90:692
gmpdata::igtp
integer(mpi), dimension(numgxy) igtp
type of XY data
Definition: mphistab.f90:675
gmpdata::lvers
integer(mpi), dimension(numgxy) lvers
version
Definition: mphistab.f90:681
gmpdata::xyplws
real(mps), dimension(10, numgxy) xyplws
additional data for GMPMS
Definition: mphistab.f90:691
gmpdata::nstr
integer(mpi), dimension(numgxy) nstr
initialization flag
Definition: mphistab.f90:674
gmpdata::nst
integer(mpi), dimension(3, numgxy) nst
counters for GMPMS
Definition: mphistab.f90:682
hmpcons
Histogram constants.
Definition: mphistab.f90:93
hmpcons::nbin
integer(mpi), parameter nbin
number of bins
Definition: mphistab.f90:98
hmpcons::numhis
integer(mpi), parameter numhis
number of histograms
Definition: mphistab.f90:97
hmpcons::nsampl
integer(mpi), parameter nsampl
number of samples for auto scaling
Definition: mphistab.f90:99
hmpdata
Histogram data.
Definition: mphistab.f90:104
hmpdata::khist
integer(mpi), dimension(numhis) khist
histgram type
Definition: mphistab.f90:111
hmpdata::xl
real(mps), dimension(6, numhis) xl
histogram binning
Definition: mphistab.f90:114
hmpdata::lun
integer(mpi) lun
unit for output
Definition: mphistab.f90:108
hmpdata::htext
character(len=60), dimension(numhis) htext
histogram text
Definition: mphistab.f90:116
hmpdata::kvers
integer(mpi), dimension(numhis) kvers
histogram version
Definition: mphistab.f90:112
hmpdata::fnhist
real(mps), dimension(nsampl, numhis) fnhist
initial data for auto scaling
Definition: mphistab.f90:113
hmpdata::dl
real(mpd), dimension(2, numhis) dl
histogram moments
Definition: mphistab.f90:115
hmpdata::jnhist
integer(mpi), dimension(5, numhis) jnhist
histogram statistics
Definition: mphistab.f90:110
hmpdata::inhist
integer(mpi), dimension(nbin, numhis) inhist
histogram (bin) data
Definition: mphistab.f90:109
mpdef
Definition of constants.
Definition: mpdef.f90:24
mpdef::mps
integer, parameter mps
single precision
Definition: mpdef.f90:37
stmamod
storage manager data.
Definition: mphistab.f90:1063
stmamod::iflc2
integer(mpi) iflc2
Definition: mphistab.f90:1070
stmamod::next
integer(mpi), dimension(:), allocatable next
Definition: mphistab.f90:1068
stmamod::tk
real(mps), dimension(:,:), allocatable tk
Definition: mphistab.f90:1067
stmamod::iflc1
integer(mpi) iflc1
Definition: mphistab.f90:1069
pfvert
subroutine pfvert(n, x)
Vertical print of floating point data.
Definition: vertpr.f90:225
pivert
subroutine pivert(n, list)
Vertical print of integer data.
Definition: vertpr.f90:185
psvert
subroutine psvert(xa, xb)
Print scale.
Definition: vertpr.f90:259
Generated on Wed Sep 18 2024 14:36:52 for Millepede-II by doxygen 1.9.3