source: trunk/readexp.tcl @ 1166

Last change on this file since 1166 was 1166, checked in by toby, 9 years ago

bring sandbox changes over to main release

  • Property svn:keywords set to Author Date Revision Id
File size: 119.9 KB
Line 
1# $Id: readexp.tcl 1166 2011-08-17 23:17:04Z toby $
2# Routines to deal with the .EXP "data structure"
3set expmap(Revision) {$Revision: 1166 $ $Date: 2011-08-17 23:17:04 +0000 (Wed, 17 Aug 2011) $}
4
5#  The GSAS data is read from an EXP file.
6#   ... reading an EXP file into an array
7# returns -1 on error
8# returns 0 if the file is old-style UNIX format (no CR/LF)
9# returns 1 if the file is 80 char/line + cr/lf
10# returns 2 if the file is sequential but not fixed-record length
11proc expload {expfile "ns {}"} {
12    # expfile is the path to the data file.
13    # ns is the namespace to place the output array (default is global)
14    if {$ns != ""} {
15        namespace eval $ns {}
16    }
17    if [catch {set fil [open "$expfile" r]}] {
18        tk_dialog .expFileErrorMsg "File Open Error" \
19                "Unable to open file $expfile" error 0 "Exit" 
20        return -1
21    }
22    fconfigure $fil -translation lf
23    set len [gets $fil line]
24    if {[string length $line] != $len} {
25        tk_dialog .expConvErrorMsg "old tcl" \
26                "You are using an old version of Tcl/Tk and your .EXP file has binary characters; run convstod or upgrade" \
27                error 0 "Exit"
28        return -1
29    }
30    catch {
31        unset ${ns}::exparray
32    }
33    if {$len > 160} {
34        set fmt 0
35        # a UNIX-type file
36        set i1 0
37        set i2 79
38        while {$i2 < $len} {
39            set nline [string range $line $i1 $i2]
40            incr i1 80
41            incr i2 80
42            set key [string range $nline 0 11]
43            set ${ns}::exparray($key) [string range $nline 12 end]
44        }
45    } else {
46        set fmt 1
47        while {$len > 0} {
48            set key [string range $line 0 11]
49            set ${ns}::exparray($key) [string range $line 12 79]
50            if {$len != 81 || [string range $line end end] != "\r"} {set fmt 2}
51            set len [gets $fil line]
52        }
53    }
54    close $fil
55    return $fmt
56}
57
58proc createexp {expfile title} {
59    global exparray expmap
60    catch {unset exparray}
61    foreach key   {"     VERSION" "      DESCR" "ZZZZZZZZZZZZ" " EXPR NPHAS"} \
62            value {"   6"         ""            "  Last EXP file record" ""} {
63        # truncate long keys & pad short ones
64        set key [string range "$key        " 0 11]
65        set exparray($key) $value
66    }
67    expinfo title set $title
68    exphistory add " created readexp.tcl [lindex $expmap(Revision) 1] [clock format [clock seconds] -format %Y-%m-%dT%T]"
69    expwrite $expfile
70}
71
72# get information out from an EXP file
73#   creates the following entries in global array expmap
74#     expmap(phaselist)     gives a list of defined phases
75#     expmap(phasetype)     gives the phase type for each defined phase
76#                           =1 nuclear; 2 mag+nuc; 3 mag; 4 macro
77#     expmap(atomlist_$p)   gives a list of defined atoms in phase $p
78#     expmap(htype_$n)      gives the GSAS histogram type for histogram (all)
79#     expmap(powderlist)    gives a list of powder histograms in use
80#     expmap(phaselist_$n)  gives a list of phases used in histogram $n
81#     expmap(nhst)          the number of GSAS histograms
82#
83proc mapexp {} {
84    global expgui expmap exparray
85    # clear out the old array
86    set expmap_Revision $expmap(Revision)
87    unset expmap
88    set expmap(Revision) $expmap_Revision
89    # apply any updates to the .EXP file
90    updateexpfile
91    # get the defined phases
92    set line [readexp " EXPR NPHAS"]
93#    if {$line == ""} {
94#       set msg "No EXPR NPHAS entry. This is an invalid .EXP file"
95#       tk_dialog .badexp "Error in EXP" $msg error 0 Exit
96#       destroy .
97#    }
98    set expmap(phaselist) {}
99    set expmap(phasetype) {}
100    # loop over phases
101    foreach iph {1 2 3 4 5 6 7 8 9} {
102        set i5s [expr {($iph - 1)*5}]
103        set i5e [expr {$i5s + 4}]
104        set flag [string trim [string range $line $i5s $i5e]]
105        if {$flag == ""} {set flag 0}
106        if $flag {
107            lappend expmap(phaselist) $iph
108            lappend expmap(phasetype) $flag
109        }
110    }
111    # get the list of defined atoms for each phase
112    foreach iph $expmap(phaselist) {
113        set expmap(atomlist_$iph) {}
114        if {[lindex $expmap(phasetype) [expr {$iph - 1}]] != 4} {
115            foreach key [array names exparray "CRS$iph  AT*A"] {
116                regexp { AT *([0-9]+)A} $key a num
117                lappend expmap(atomlist_$iph) $num
118            }
119        } else {
120            foreach key [array names exparray "CRS$iph  AT*"] {
121                scan [string range $key 8 11] %x atm
122                lappend expmap(atomlist_$iph) $atm
123            }
124        }
125        # note that sometimes an .EXP file contains more atoms than are actually defined
126        # drop the extra ones
127        set expmap(atomlist_$iph) [lsort -integer $expmap(atomlist_$iph)]
128        set natom [phaseinfo $iph natoms]
129        if {$natom != [llength $expmap(atomlist_$iph)]} {
130            set expmap(atomlist_$iph) [lrange $expmap(atomlist_$iph) 0 [expr {$natom-1}]]
131        }
132    }
133    # now get the histogram types
134    set expmap(nhst) [string trim [readexp { EXPR  NHST }]]
135    set n 0
136    set expmap(powderlist) {}
137    for {set i 0} {$i < $expmap(nhst)} {incr i} {
138        set ihist [expr {$i + 1}]
139        if {[expr {$i % 12}] == 0} {
140            incr n
141            set line [readexp " EXPR  HTYP$n"]
142            if {$line == ""} {
143                set msg "No HTYP$n entry for Histogram $ihist. This is an invalid .EXP file"
144                tk_dialog .badexp "Error in readexp" $msg error 0 Exit
145            }
146            set j 0
147        } else {
148            incr j
149        }
150        set expmap(htype_$ihist) [string range $line [expr 2+5*$j] [expr 5*($j+1)]]
151        # is this a dummy histogram?
152        if {$ihist <=9} {
153            set key "HST  $ihist DUMMY"
154        } else {
155            set key "HST $ihist DUMMY"
156        }
157        # at least for now, ignore non-powder histograms
158        if {[string range $expmap(htype_$ihist) 0 0] == "P" && \
159                [string range $expmap(htype_$ihist) 3 3] != "*"} {
160            if {[existsexp $key]} {
161                set expmap(htype_$ihist) \
162                        [string range $expmap(htype_$ihist) 0 2]D
163            }
164            lappend expmap(powderlist) $ihist
165        }
166    }
167
168    # now process powder histograms
169    foreach ihist $expmap(powderlist) {
170        # make a 2 digit key -- hh
171        if {$ihist < 10} {
172            set hh " $ihist"
173        } else {
174            set hh $ihist
175        }
176        set line [readexp "HST $hh NPHAS"]
177        if {$line == ""} {
178            set msg "No NPHAS entry for Histogram $ihist. This is an invalid .EXP file"
179            tk_dialog .badexp "Error in readexp" $msg error 0 Exit
180        }
181        set expmap(phaselist_$ihist) {}
182        # loop over phases
183        foreach iph {1 2 3 4 5 6 7 8 9} {
184            set i5s [expr {($iph - 1)*5}]
185            set i5e [expr {$i5s + 4}]
186            set flag [string trim [string range $line $i5s $i5e]]
187            if {$flag == ""} {set flag 0}
188            if $flag {lappend expmap(phaselist_$ihist) $iph}
189        }
190    }
191    # load the constrained parameters
192    atom_constraint_load
193    set expgui(mapstat) 1
194}
195
196# this routine is called to update changes in the .EXP file
197proc updateexpfile {} {
198    global exparray
199    # change "CRSx  ODFxxx" records to "CRSx  OD xxx" records
200    # needed by the June 8, 2005 GSAS release
201    set ODFlist [array names exparray "CRS?  ODF*"]
202    if {[llength $ODFlist] > 0} {
203        catch {incr ::expgui(changed)}
204        foreach key $ODFlist {
205            regsub "ODF" $key "OD " newkey
206            set exparray($newkey) $exparray($key)
207            unset exparray($key)
208        }
209    }
210}
211
212
213# return the value for a ISAM key
214proc readexp {key} {
215    global exparray
216    # truncate long keys & pad short ones
217    set key [string range "$key        " 0 11]
218    if [catch {set val $exparray($key)}] {
219        global expgui
220        if $expgui(debug) {puts "Error accessing record $key"}
221        return ""
222    }
223    return $val
224}
225
226# return the number of records matching ISAM key (may contain wildcards)
227proc existsexp {key} {
228    global exparray
229    # key can contain wild cards so don't pad
230    return [llength [array names exparray  $key]]
231}
232
233
234# replace a section of the exparray with $value
235#   replace $char characters starting at character $start (numbered from 1)
236proc setexp {key value start chars} {
237    global exparray
238    # truncate long keys & pad short ones
239    set key [string range "$key        " 0 11]
240    if [catch {set exparray($key)}] {
241        global expgui
242        if $expgui(debug) {puts "Error accessing record $key"}
243        return ""
244    }
245
246    # pad value to $chars
247    set l0 [expr {$chars - 1}]
248    set value [string range "$value                                           " 0 $l0]
249
250    if {$start == 1} {
251        set ret {}
252        set l1 $chars
253    } else {
254        set l0 [expr {$start - 2}]
255        set l1 [expr {$start + $chars - 1}]
256        set ret [string range $exparray($key) 0 $l0]
257    }
258    append ret $value [string range $exparray($key) $l1 end]
259    set exparray($key) $ret
260}
261
262proc makeexprec {key} {
263    global exparray
264    # truncate long keys & pad short ones
265    set key [string range "$key        " 0 11]
266    if [catch {set exparray($key)}] {
267        # set to 68 blanks
268        set exparray($key) [format %68s " "]
269    }
270}
271
272# delete an exp record
273# returns 1 if OK; 0 if not found
274proc delexp {key} {
275    global exparray
276    # truncate long keys & pad short ones
277    set key [string range "$key        " 0 11]
278    if [catch {unset exparray($key)}] {
279        return 0
280    }
281    return 1
282}
283
284# test an argument if it is a valid number; reform the number to fit
285proc validreal {val length decimal} {
286    upvar $val value
287    # is this a number?
288    if [catch {expr {$value}}] {return 0}
289    if [catch {
290        # how many digits are needed to the left of the decimal?
291        set sign 0
292        if {$value > 0} {
293            set digits [expr {1 + int(log10($value))}]
294        } elseif {$value < 0} {
295            set digits [expr {1 + int(log10(-$value))}]
296            set sign 1
297        } else {
298            set digits 1
299        }
300        if {$digits <= 0} {set digits 1}
301        if {$digits + $sign >= $length} {
302            # the number is much too big -- use exponential notation
303            set decimal [expr {$length - 6 - $sign}]
304            # drop more decimal places, as needed
305            set tmp [format "%${length}.${decimal}E" $value]
306            while {[string length $tmp] > $length && $decimal >= 0} {
307                incr decimal -1
308                set tmp [format "%${length}.${decimal}E" $value]
309            }
310        } elseif {$digits + $sign >= $length - $decimal} {
311            # we will have to trim the number of decimal digits
312            set decimal [expr {$length - $digits - $sign - 1}]
313            set tmp [format "%#.${decimal}f" $value]
314        } elseif {abs($value) < pow(10,2-$decimal) && $length > 6 && $value != 0} {
315            # for small values, switch to exponential notation (2-$decimal -> three sig figs.)
316            set decimal [expr {$length - 6 - $sign}]
317            # drop more decimal places, as needed
318            set tmp [format "%${length}.${decimal}E" $value]
319            while {[string length $tmp] > $length && $decimal >= 0} {
320                incr decimal -1
321                set tmp [format "%${length}.${decimal}E" $value]
322            }
323        } else {
324            # use format as specified
325            set tmp [format "%${length}.${decimal}f" $value]
326        }
327        set value $tmp
328    } errmsg] {return 0}
329    return 1
330}
331
332# test an argument if it is a valid integer; reform the number into
333# an integer, if appropriate -- be sure to pass the name of the variable not the value
334proc validint {val length} {
335    upvar $val value
336    # FORTRAN type assumption: blank is 0
337    if {$value == ""} {set value 0}
338    if [catch {
339        set tmp [expr {round($value)}]
340        if {$tmp != $value} {return 0}
341        set value [format "%${length}d" $tmp]
342    }] {return 0}
343    return 1
344}
345
346# process history information
347#    action == last
348#       returns number and value of last record
349#    action == add
350#
351proc exphistory {action "value 0"} {
352    global exparray
353    if {$action == "last"} {
354        set key [lindex [lsort -decreasing [array names exparray *HSTRY*]] 0]
355        if {$key == ""} {return ""}
356        return [list [string trim [string range $key 9 end]] $exparray($key)]
357    } elseif {$action == "add"} {
358        set key [lindex [lsort -decreasing [array names exparray *HSTRY*]] 0]
359        if {$key == ""} {
360            set index 1
361        } else {
362            set index [string trim [string range $key 9 end]]
363            if {$index != "***"} {
364                if {$index < 999} {incr index}
365                set key [format "    HSTRY%3d" $index]
366                set exparray($key) $value
367            }
368        }
369        set key [format "    HSTRY%3d" $index]
370        set exparray($key) $value
371    }
372}
373# get overall info
374#   parm:
375#     print     -- GENLES print option (*)
376#     cycles    -- number of GENLES cycles (*)
377#     title     -- the overall title (*)
378#     convg     -- convergence criterion: -200 to 200 (*)
379#     marq      -- Marquardt damping factor: 1.0 to 9.99 (*)
380#     mbw       -- LS matrix bandwidth; =0 for full matrix (*)
381proc expinfo {parm "action get" "value {}"} {
382    switch ${parm}-$action {
383        title-get {
384            return [string trim [readexp "      DESCR"]]
385        }
386        title-set {
387            setexp "      DESCR" "  $value" 2 68
388        }
389        cycles-get {
390            return [string trim [cdatget MXCY]]
391        }
392        cycles-set {
393            if ![validint value 1] {return 0}
394            cdatset MXCY [format %4d $value]
395        }
396        cyclesrun-get {
397            set cycle -1
398            regexp {.*cycles run *([0-9]*) } [readexp "  GNLS  RUN"] x cycle
399            return $cycle
400        }
401        print-get {
402            set print [string trim [cdatget PRNT]]
403            if {$print != ""} {return $print}
404            return 0
405        }
406        print-set {
407            if ![validint value 1] {return 0}
408            cdatset PRNT [format %4d $value]
409        }
410        convg-get {
411            set cvg [string trim [cdatget CVRG]]
412            if {$cvg == ""} {return -200}
413            if [catch {expr {$cvg}}] {return -200}
414            return $cvg
415        }
416        convg-set {
417            if ![validint value 1] {return 0}
418            set value [expr {-200>$value?-200:$value}]
419            set value [expr {200<$value?200:$value}]
420            cdatset CVRG [format %4d $value]
421        }
422        marq-get {
423            set mar [string trim [cdatget MARQ]]
424            if {$mar == ""} {return 1.0}
425            if [catch {expr $mar}] {return 1.}
426            return $mar
427        }
428        marq-set {
429            if [catch {
430                set value [expr {1.0>$value?1.0:$value}]
431                set value [expr {9.99<$value?9.99:$value}]
432            }] {return 0}
433            if ![validreal value 4 2] {return 0}
434            cdatset MARQ $value
435        }
436        mbw-get {
437            set mbw [string trim [cdatget MBW]]
438            if {$mbw == ""} {return 0}
439            if [catch {expr $mbw}] {return 0}
440            return $mbw
441        }
442        mbw-set {
443            if ![validint value 1] {return 0}
444            if {$value < 0} {return 0}
445            cdatset MBW [format %5d $value]
446        }
447        default {
448            set msg "Unsupported expinfo access: parm=$parm action=$action"
449            tk_dialog .badexp "Error in readexp" $msg error 0 Exit
450        }
451    }
452    return 1
453}
454
455proc cdatget {key} {
456    foreach i {1 2 3 4 5 6 7 8 9} {
457        if {[existsexp "  GNLS CDAT$i"] == 0} break
458        set line [readexp "  GNLS CDAT$i"]
459        if {$line == {}} break
460        foreach i1 {2 10 18 26 34 42 50 58 66} \
461                i2 {9 17 25 33 41 49 57 65 73} {
462            set item [string range $line $i1 $i2]
463            if {[string trim $item] == {}} continue
464            if [regexp "${key}(.*)" $item a b] {return $b}
465        }
466    }
467    return {}
468}
469
470proc cdatset {key value} {
471    # round 1 see if we can find the string
472    foreach i {1 2 3 4 5 6 7 8 9} {
473        set line [readexp "  GNLS CDAT$i"]
474        if {$line == {}} break
475        foreach i1 {2 10 18 26 34 42 50 58 66} \
476                i2 {9 17 25 33 41 49 57 65 73} {
477            set item [string range $line $i1 $i2]
478            if {[string trim $item] == {}} continue
479            if [regexp "${key}(.*)" $item a b] {
480                # found it now replace it
481                incr i1
482                setexp "  GNLS CDAT$i" "${key}${value}" $i1 8
483                return
484            }
485        }
486    }
487    # not found, take the 1st blank space, creating a card if needed
488    foreach i {1 2 3 4 5 6 7 8 9} {
489        set line [readexp "  GNLS CDAT$i"]
490        if {$line == {}} {makeexprec "  GNLS CDAT$i"}
491        foreach i1 {2 10 18 26 34 42 50 58 66} \
492                i2 {9 17 25 33 41 49 57 65 73} {
493            set item [string range $line $i1 $i2]
494            if {[string trim $item] == {}} {
495                # found a blank space: now replace it
496                incr i1
497                setexp "  GNLS CDAT$i" "${key}${value}" $i1 8
498                return
499            }
500        }
501    }
502    return {}
503}
504
505proc disagldat_get {phase} {
506    set key "  DSGL CDAT$phase"
507    if {[existsexp $key] == 0} {return "{none} {none}"}
508    set line [readexp $key]
509    set i1 2
510    # read atom-atom distance parameter
511    set dist {}
512    set item [string range $line $i1 [expr {$i1+3}]]
513    if {$item == "DMAX"} {
514        set val [string range $line [expr {$i1+4}] [expr {$i1+11}]]
515        set dist [string trim $val]
516        incr i1 13
517    } else {
518        set dist "radii"
519        incr i1 5
520    }
521    # read atom-atom-atom angle parameter
522    set ang {}
523    set item [string range $line $i1 [expr {$i1+3}]]
524    if {$item == "DAGL"} {
525        set val [string range $line [expr {$i1+4}] [expr {$i1+11}]]
526        set ang [string trim $val]
527        incr i1 13
528    } else {
529        set ang "radii"
530        incr i1 5
531    }
532    # note there are two more parameters, NOFO/FORA & ONCR/DFLT, but they are not being processed yet
533    return "$dist $ang"
534}
535
536# get phase information: phaseinfo phase parm action value
537#   phase: 1 to 9 (as defined)
538#   parm:
539#     name -- phase name
540#     natoms -- number of atoms (*)
541#     a b c alpha beta gamma -- cell parameters (*)
542#     cellref -- refinement flag for the unit cell(*)
543#     celldamp  -- damping for the unit cell refinement (*)
544#     spacegroup -- space group symbol (*)
545#     ODForder -- spherical harmonic order (*)
546#     ODFsym   -- sample symmetry (0-3) (*)
547#     ODFdampA -- damping for angles (*)
548#     ODFdampC -- damping for coefficients (*)
549#     ODFomega -- omega oriention angle (*)
550#     ODFchi -- chi oriention angle (*)
551#     ODFphi -- phi oriention angle (*)
552#     ODFomegaRef -- refinement flag for omega (*)
553#     ODFchiRef -- refinement flag for chi (*)
554#     ODFphiRef -- refinement flag for phi (*)
555#     ODFterms -- a list of the {l m n} values for each ODF term (*)
556#     ODFcoefXXX -- the ODF coefficient for for ODF term XXX (*)
557#     ODFRefcoef -- refinement flag for ODF terms (*)
558#     DistCalc   -- returns "radii", "none" or a number (*)
559#                   none: no distance or angle computation for the phase
560#                   radii: computation will be done by sums of radii
561#                          (see AtmTypInfo and DefAtmTypInfo)
562#                   other: a distance specifing the maximum distance
563#     AngCalc    -- returns "radii", "none" or a number (*)
564#                   none: no distance or angle computation for the phase
565#                   radii: computation will be done by sums of radii
566#                          (see AtmTypInfo and DefAtmTypInfo)
567#                   other: a distance specifing the maximum distance
568#  action: get (default) or set
569#  value: used only with set
570#  * =>  read+write supported
571proc phaseinfo {phase parm "action get" "value {}"} {
572    switch -glob ${parm}-$action {
573
574        name-get {
575            return [string trim [readexp "CRS$phase    PNAM"]]
576        }
577
578        name-set {
579            setexp "CRS$phase    PNAM" " $value" 2 68
580        }
581
582        spacegroup-get {
583            return [string trim [readexp "CRS$phase  SG SYM"]]
584        }
585
586        spacegroup-set {
587            setexp "CRS$phase  SG SYM" " $value" 2 68
588        }
589
590        natoms-get {
591            return [string trim [readexp "CRS$phase   NATOM"]]     
592        }
593
594        natoms-set {
595            if ![validint value 5] {return 0}
596            setexp "CRS$phase   NATOM" $value 1 5
597        }
598
599        a-get {
600           return [string trim [string range [readexp "CRS$phase  ABC"] 0 9]]
601        }
602        b-get {
603           return [string trim [string range [readexp "CRS$phase  ABC"] 10 19]]
604        }
605        c-get {
606           return [string trim [string range [readexp "CRS$phase  ABC"] 20 29]]
607        }
608        alpha-get {
609           return [string trim [string range [readexp "CRS$phase  ANGLES"] 0 9]]
610        }
611        beta-get {
612           return [string trim [string range [readexp "CRS$phase  ANGLES"] 10 19]]
613        }
614        gamma-get {
615           return [string trim [string range [readexp "CRS$phase  ANGLES"] 20 29]]
616        }
617
618        a-set {
619            if ![validreal value 10 6] {return 0}
620            setexp "CRS$phase  ABC" $value 1 10             
621        }
622        b-set {
623            if ![validreal value 10 6] {return 0}
624            setexp "CRS$phase  ABC" $value 11 10           
625        }
626        c-set {
627            if ![validreal value 10 6] {return 0}
628            setexp "CRS$phase  ABC" $value 21 10           
629        }
630        alpha-set {
631            if ![validreal value 10 4] {return 0}
632            setexp "CRS$phase  ANGLES" $value 1 10         
633        }
634        beta-set {
635            if ![validreal value 10 4] {return 0}
636            setexp "CRS$phase  ANGLES" $value 11 10         
637        }
638        gamma-set {
639            if ![validreal value 10 4] {return 0}
640            setexp "CRS$phase  ANGLES" $value 21 10         
641        }
642        cellref-get {
643            if {[string toupper [string range [readexp "CRS$phase  ABC"] 34 34]] == "Y"} {
644                return 1
645            }
646            return 0
647        }
648        cellref-set {
649            if $value {
650                setexp "CRS$phase  ABC" "Y" 35 1
651            } else {
652                setexp "CRS$phase  ABC" "N" 35 1
653            }       
654        }
655        celldamp-get {
656            set val [string range [readexp "CRS$phase  ABC"] 39 39]
657            if {$val == " "} {return 0}
658            return $val
659        }
660        celldamp-set {
661            setexp "CRS$phase  ABC" $value 40 1
662        }
663
664        ODForder-get {
665            set val [string trim [string range [readexp "CRS$phase  OD "] 0 4]]
666            if {$val == " "} {return 0}
667            return $val
668        }
669        ODForder-set {
670            if ![validint value 5] {return 0}
671            setexp "CRS$phase  OD " $value 1 5
672        }
673        ODFsym-get {
674            set val [string trim [string range [readexp "CRS$phase  OD "] 10 14]]
675            if {$val == " "} {return 0}
676            return $val
677        }
678        ODFsym-set {
679            if ![validint value 5] {return 0}
680            setexp "CRS$phase  OD " $value 11 5
681        }
682        ODFdampA-get {
683            set val [string range [readexp "CRS$phase  OD "] 24 24]
684            if {$val == " "} {return 0}
685            return $val
686        }
687        ODFdampA-set {
688            setexp "CRS$phase  OD " $value 25 1
689        }
690        ODFdampC-get {
691            set val [string range [readexp "CRS$phase  OD "] 29 29]
692            if {$val == " "} {return 0}
693            return $val
694        }
695        ODFdampC-set {
696            setexp "CRS$phase  OD " $value 30 1
697        }
698        ODFomegaRef-get {
699            if {[string toupper [string range [readexp "CRS$phase  OD "] 16 16]] == "Y"} {
700                return 1
701            }
702            return 0
703        }
704        ODFomegaRef-set {
705            if $value {
706                setexp "CRS$phase  OD " "Y" 17 1
707            } else {
708                setexp "CRS$phase  OD " "N" 17 1
709            }       
710        }
711        ODFchiRef-get {
712            if {[string toupper [string range [readexp "CRS$phase  OD "] 17 17]] == "Y"} {
713                return 1
714            }
715            return 0
716        }
717        ODFchiRef-set {
718            if $value {
719                setexp "CRS$phase  OD " "Y" 18 1
720            } else {
721                setexp "CRS$phase  OD " "N" 18 1
722            }       
723        }
724        ODFphiRef-get {
725            if {[string toupper [string range [readexp "CRS$phase  OD "] 18 18]] == "Y"} {
726                return 1
727            }
728            return 0
729        }
730        ODFphiRef-set {
731            if $value {
732                setexp "CRS$phase  OD " "Y" 19 1
733            } else {
734                setexp "CRS$phase  OD " "N" 19 1
735            }       
736        }
737        ODFcoef*-get {
738            regsub ODFcoef $parm {} term
739            set k [expr {($term+5)/6}]
740            if {$k <= 9} {
741                set k "  $k"
742            } elseif {$k <= 99} {
743                set k " $k"
744            }
745            set j [expr {(($term-1) % 6)+1}]
746            set lineB [readexp "CRS$phase  OD${k}B"]
747            set j0 [expr { ($j-1) *10}]
748            set j1 [expr {$j0 + 9}]
749            set val [string trim [string range $lineB $j0 $j1]]
750            if {$val == ""} {return 0.0}
751            return $val
752        }
753        ODFcoef*-set {
754            regsub ODFcoef $parm {} term
755            if ![validreal value 10 3] {return 0}
756            set k [expr {($term+5)/6}]
757            if {$k <= 9} {
758                set k "  $k"
759            } elseif {$k <= 99} {
760                set k " $k"
761            }
762            set j [expr {(($term-1) % 6)+1}]
763            set col [expr { ($j-1)*10 + 1}]
764            setexp "CRS$phase  OD${k}B" $value $col 10
765        }
766        ODFRefcoef-get {
767            if {[string toupper [string range [readexp "CRS$phase  OD "] 19 19]] == "Y"} {
768                return 1
769            }
770            return 0
771        }
772        ODFRefcoef-set {
773            if $value {
774                setexp "CRS$phase  OD " "Y" 20 1
775            } else {
776                setexp "CRS$phase  OD " "N" 20 1
777            }       
778        }
779        ODFomega-get {
780           return [string trim [string range [readexp "CRS$phase  OD "] 30 39]]
781        }
782        ODFchi-get {
783           return [string trim [string range [readexp "CRS$phase  OD "] 40 49]]
784        }
785        ODFphi-get {
786           return [string trim [string range [readexp "CRS$phase  OD "] 50 59]]
787        }
788        ODFomega-set {
789            if ![validreal value 10 4] {return 0}
790            setexp "CRS$phase  OD " $value 31 10
791        }
792        ODFchi-set {
793            if ![validreal value 10 4] {return 0}
794            setexp "CRS$phase  OD " $value 41 10
795        }
796        ODFphi-set {
797            if ![validreal value 10 4] {return 0}
798            setexp "CRS$phase  OD " $value 51 10
799        }
800
801        ODFterms-get {
802            set vallist {}
803            set val [string trim [string range [readexp "CRS$phase  OD "] 5 9]]
804            for {set i 1} {$i <= $val} {incr i 6} {
805                set k [expr {1+($i-1)/6}]
806                if {$k <= 9} {
807                    set k "  $k"
808                } elseif {$k <= 99} {
809                    set k " $k"
810                }
811                set lineA [readexp "CRS$phase  OD${k}A"]
812                set k 0
813                for {set j $i} {$j <= $val && $j < $i+6} {incr j} {
814                    set j0 [expr {($k)*10}]
815                    set j1 [expr {$j0 + 9}]
816                    lappend vallist [string trim [string range $lineA $j0 $j1]]
817                    incr k
818                }
819            }
820            return $vallist
821        }
822        ODFterms-set {
823            set key "CRS$phase  OD    "
824            if {![existsexp $key]} {
825                makeexprec $key
826                set oldlen 0
827            } else {
828                set oldlen [string trim [string range [readexp $key] 5 9]]
829            }
830            set len [llength $value]
831            if ![validint len 5] {return 0}
832            setexp $key $len 6 5
833            set j 0
834            set k 0
835            foreach item $value {
836                incr j
837                if {$j % 6 == 1} {
838                    incr k
839                    if {$k <= 9} {
840                        set k "  $k"
841                    } elseif {$k <= 99} {
842                        set k " $k"
843                    }
844                    set col 1
845                    set keyA "CRS$phase  OD${k}A"
846                    set keyB "CRS$phase  OD${k}B"
847                    if {![existsexp $keyA]} {
848                        makeexprec $keyA
849                        makeexprec $keyB
850                    }
851                }
852                set col1 [expr {$col + 1}]
853                foreach n [lrange $item 0 2] {
854                    if ![validint n 3] {return 0}
855                    setexp $keyA $n $col1 3
856                    incr col1 3
857                }
858                incr col 10
859            }
860            for {incr j} {$j <= $oldlen} {incr j} {
861                if {$j % 6 == 1} {
862                    incr k
863                    if {$k <= 9} {
864                        set k "  $k"
865                    } elseif {$k <= 99} {
866                        set k " $k"
867                    }
868                    set col 1
869                    set keyA "CRS$phase  OD${k}A"
870                    set keyB "CRS$phase  OD${k}B"
871                    delexp $keyA
872                    delexp $keyB
873                }
874                if {[existsexp $keyA]} {
875                    setexp $keyA "          " $col 10
876                    setexp $keyB "          " $col 10
877                }
878                incr col 10
879            }
880        }
881        DistCalc-get {
882            set val [disagldat_get $phase]
883            return [lindex $val 0]
884        }
885        DistCalc-set {
886            set key "  DSGL CDAT$phase"
887            # for none delete the record & thats all folks
888            if {$value == "none"} {
889                catch {unset ::exparray($key)}
890                return
891            }
892            if {[existsexp $key] == 0} {
893                makeexprec $key
894            }
895            set line [readexp $key]
896            if {[string trim $line] == ""} {
897                # blank set to defaults
898                set line [string replace $line 2 15 "DRAD ARAD NOFO"]
899            }
900            if {$value == "radii"} {
901                if {[string range $line 2 5] == "DMAX"} {
902                    set line [string replace $line 2 13 "DRAD"]
903                } else {
904                    set line [string replace $line 2 5 "DRAD"]
905                }
906            } else {
907                if ![validreal value 8 2] {return 0}
908                if {[string range $line 2 5] == "DMAX"} {
909                    set line [string replace $line 6 13 $value]
910                } else {
911                    set line [string replace $line 2 5 "DMAX"]
912                    set line [string replace $line 6 6 "$value "]
913                }
914            }
915            setexp $key $line 0 68
916        }
917        AngCalc-get {
918            set val [disagldat_get $phase]
919            return [lindex $val 1]
920        }
921        AngCalc-set {
922            set key "  DSGL CDAT$phase"
923            # for none delete the record & thats all folks
924            if {$value == "none"} {
925                catch {unset ::exparray($key)}
926                return
927            }
928            if {[existsexp $key] == 0} {
929                makeexprec $key
930            }
931            set line [readexp $key]
932            if {[string trim $line] == ""} {
933                # blank set to defaults
934                set line [string replace $line 2 15 "DRAD ARAD NOFO"]
935            }
936            if {[string range $line 2 5] == "DMAX"} {
937                set i2 8
938            } else {
939                set i2 0
940            }
941            if {$value == "radii"} {
942                if {[string range $line [expr {$i2+7}] [expr {$i2+10}]] == "DAGL"} {
943                    set line [string replace $line [expr {$i2+7}] [expr {$i2+18}] "ARAD"]
944                } else {
945                    set line [string replace $line [expr {$i2+7}] [expr {$i2+10}] "ARAD"]
946                }
947            } else {
948                if ![validreal value 8 2] {return 0}
949                if {[string range $line [expr {$i2+7}] [expr {$i2+10}]] == "DAGL"} {
950                    set line [string replace $line [expr {$i2+11}] [expr {$i2+18}] $value]
951                } else {
952                    set line [string replace $line [expr {$i2+7}] [expr {$i2+10}] "DAGL"]
953                    set line [string replace $line [expr {$i2+11}] [expr {$i2+11}] "$value "]
954                }
955            }
956            setexp $key $line 0 68
957        }
958        default {
959            set msg "Unsupported phaseinfo access: parm=$parm action=$action"
960            tk_dialog .badexp "Error in readexp" $msg error 0 Exit
961        }
962    }
963    return 1
964}
965
966
967
968# get atom information: atominfo phase atom parm action value
969#   phase: 1 to 9 (as defined)
970#   atom: a valid atom number [see expmap(atomlist_$phase)]
971#      Note that atom and phase can be paired lists, but if there are extra
972#      entries in the atoms list, the last phase will be repeated.
973#      so that atominfo 1 {1 2 3} xset 1
974#               will set the xflag for atoms 1-3 in phase 1
975#      but atominfo {1 2 3} {1 1 1} xset 1
976#               will set the xflag for atoms 1 in phase 1-3
977#   parm:
978#     type -- element code
979#     mult -- atom multiplicity
980#     label -- atom label (*)
981#     x y z -- coordinates (*)
982#     frac --  occupancy (*)
983#     temptype -- I or A for Isotropic/Anisotropic (*)
984#     Uiso  -- Isotropic temperature factor (*)
985#     U11  -- Anisotropic temperature factor (*)
986#     U22  -- Anisotropic temperature factor (*)
987#     U33  -- Anisotropic temperature factor (*)
988#     U12  -- Anisotropic temperature factor (*)
989#     U13  -- Anisotropic temperature factor (*)
990#     U23  -- Anisotropic temperature factor (*)
991#     xref/xdamp -- refinement flag/damping value for the coordinates (*)
992#     uref/udamp -- refinement flag/damping value for the temperature factor(s)  (*)
993#     fref/fdamp -- refinement flag/damping value for the occupancy (*)
994#  action: get (default) or set
995#  value: used only with set
996#  * =>  read+write supported
997proc atominfo {phaselist atomlist parm "action get" "value {}"} {
998    foreach phase $phaselist atom $atomlist {
999        if {$phase == ""} {set phase [lindex $phaselist end]}
1000        if {$atom < 10} {
1001            set key "CRS$phase  AT  $atom"
1002        } elseif {$atom < 100} {
1003            set key "CRS$phase  AT $atom"
1004        } else {
1005            set key "CRS$phase  AT$atom"
1006        }
1007        switch -glob ${parm}-$action {
1008            type-get {
1009                return [string trim [string range [readexp ${key}A] 2 9] ]
1010            }
1011            mult-get {
1012                return [string trim [string range [readexp ${key}A] 58 61] ]
1013            }
1014            label-get {
1015                return [string trim [string range [readexp ${key}A] 50 57] ]
1016            }
1017            label-set {
1018                setexp ${key}A $value 51 8
1019            }
1020            temptype-get {
1021                return [string trim [string range [readexp ${key}B] 62 62] ]
1022            }
1023            temptype-set {
1024                if {$value == "A"} {
1025                    setexp ${key}B A 63 1
1026                    # copy the Uiso to the diagonal terms
1027                    set Uiso [string range [readexp ${key}B] 0 9]
1028                    foreach value [CalcAniso $phase $Uiso] \
1029                            col {1 11 21 31 41 51} {
1030                        validreal value 10 6
1031                        setexp ${key}B $value $col 10
1032                    }
1033                } else {
1034                    setexp ${key}B I 63 1
1035                    set value 0.0
1036                    catch {
1037                        # get the trace
1038                        set value [expr {( \
1039                                [string range [readexp ${key}B] 0 9] + \
1040                                [string range [readexp ${key}B] 10 19] + \
1041                                [string range [readexp ${key}B] 20 29])/3.}]
1042                    }
1043                    validreal value 10 6
1044                    setexp ${key}B $value 1 10
1045                    # blank out the remaining terms
1046                    set value " "
1047                    setexp ${key}B $value 11 10
1048                    setexp ${key}B $value 21 10
1049                    setexp ${key}B $value 31 10
1050                    setexp ${key}B $value 41 10
1051                    setexp ${key}B $value 51 10
1052                }
1053            }
1054            x-get {
1055                return [string trim [string range [readexp ${key}A] 10 19] ]
1056            }
1057            x-set {
1058                if ![validreal value 10 6] {return 0}
1059                setexp ${key}A $value 11 10
1060            }
1061            y-get {
1062                return [string trim [string range [readexp ${key}A] 20 29] ]
1063            }
1064            y-set {
1065                if ![validreal value 10 6] {return 0}
1066                setexp ${key}A $value 21 10
1067            }
1068            z-get {
1069                return [string trim [string range [readexp ${key}A] 30 39] ]
1070            }
1071            z-set {
1072                if ![validreal value 10 6] {return 0}
1073                setexp ${key}A $value 31 10
1074            }
1075            frac-get {
1076                return [string trim [string range [readexp ${key}A] 40 49] ]
1077            }
1078            frac-set {
1079                if ![validreal value 10 6] {return 0}
1080                setexp ${key}A $value 41 10
1081            }
1082            U*-get {
1083                regsub U $parm {} type
1084                if {$type == "iso" || $type == "11"} {
1085                    return [string trim [string range [readexp ${key}B] 0 9] ]
1086                } elseif {$type == "22"} {
1087                    return [string trim [string range [readexp ${key}B] 10 19] ]
1088                } elseif {$type == "33"} {
1089                    return [string trim [string range [readexp ${key}B] 20 29] ]
1090                } elseif {$type == "12"} {
1091                    return [string trim [string range [readexp ${key}B] 30 39] ]
1092                } elseif {$type == "13"} {
1093                    return [string trim [string range [readexp ${key}B] 40 49] ]
1094                } elseif {$type == "23"} {
1095                    return [string trim [string range [readexp ${key}B] 50 59] ]
1096                }
1097            }
1098            U*-set {
1099                if ![validreal value 10 6] {return 0}
1100                regsub U $parm {} type
1101                if {$type == "iso" || $type == "11"} {
1102                    setexp ${key}B $value 1 10
1103                } elseif {$type == "22"} {
1104                    setexp ${key}B $value 11 10
1105                } elseif {$type == "33"} {
1106                    setexp ${key}B $value 21 10
1107                } elseif {$type == "12"} {
1108                    setexp ${key}B $value 31 10
1109                } elseif {$type == "13"} {
1110                    setexp ${key}B $value 41 10
1111                } elseif {$type == "23"} {
1112                    setexp ${key}B $value 51 10
1113                }
1114            }
1115            xref-get {
1116                if {[string toupper [string range [readexp ${key}B] 64 64]] == "X"} {
1117                    return 1
1118                }
1119                return 0
1120            }
1121            xref-set {
1122                if $value {
1123                    setexp ${key}B "X" 65 1
1124                } else {
1125                    setexp ${key}B " " 65 1
1126                }           
1127            }
1128            xdamp-get {
1129                set val [string range [readexp ${key}A] 64 64]
1130                if {$val == " "} {return 0}
1131                return $val
1132            }
1133            xdamp-set {
1134                setexp ${key}A $value 65 1
1135            }
1136            fref-get {
1137                if {[string toupper [string range [readexp ${key}B] 63 63]] == "F"} {
1138                    return 1
1139                }
1140                return 0
1141            }
1142            fref-set {
1143                if $value {
1144                    setexp ${key}B "F" 64 1
1145                } else {
1146                    setexp ${key}B " " 64 1
1147                }           
1148            }
1149            fdamp-get {
1150                set val [string range [readexp ${key}A] 63 63]
1151                if {$val == " "} {return 0}
1152                return $val
1153            }
1154            fdamp-set {
1155                setexp ${key}A $value 64 1
1156            }
1157
1158            uref-get {
1159                if {[string toupper [string range [readexp ${key}B] 65 65]] == "U"} {
1160                    return 1
1161                }
1162                return 0
1163            }
1164            uref-set {
1165                if $value {
1166                    setexp ${key}B "U" 66 1
1167                } else {
1168                    setexp ${key}B " " 66 1
1169                }           
1170            }
1171            udamp-get {
1172                set val [string range [readexp ${key}A] 65 65]
1173                if {$val == " "} {return 0}
1174                return $val
1175            }
1176            udamp-set {
1177                setexp ${key}A $value 66 1
1178            }
1179            default {
1180                set msg "Unsupported atominfo access: parm=$parm action=$action"
1181                tk_dialog .badexp "Error in readexp" $msg error 0 Exit
1182            }
1183        }
1184    }
1185    return 1
1186}
1187
1188# get macromolecular atom information: mmatominfo phase atom parm action value
1189#   phase: 1 (at present only one mm phase can be defined)
1190#   atom: a valid atom number [see expmap(atomlist_$phase)]
1191#      Note that atoms can be lists
1192#      so that mmatominfo 1 {1 2 3} xset 1
1193#               will set the xflag for atoms 1-3 in phase 1
1194#   parm:
1195#     type -- element code
1196#     frac --  occupancy (*)
1197#     x y z -- coordinates (*)
1198#     Uiso  -- Isotropic temperature factor (*)
1199#     label -- atom label (*)
1200#     residue -- residue label (*)
1201#     group -- group label (*)
1202#     resnum -- residue number (*)
1203#     xref/xdamp -- refinement flag/damping value for the coordinates (*)
1204#     uref/udamp -- refinement flag/damping value for the temperature factor(s)  (*)
1205#     fref/fdamp -- refinement flag/damping value for the occupancy (*)
1206#  action: get (default) or set
1207#  value: used only with set
1208#  * =>  read+write supported
1209proc mmatominfo {phaselist atomlist parm "action get" "value {}"} {
1210    foreach phase $phaselist atom $atomlist {
1211        if {$phase == ""} {set phase [lindex $phaselist end]}
1212        set num [string toupper [format %.4x $atom]]
1213        set key "CRS$phase  AT$num"
1214        switch -glob ${parm}-$action {
1215            type-get {
1216                return [string trim [string range [readexp ${key}] 2 9] ]
1217            }
1218            frac-get {
1219                return [string trim [string range [readexp ${key}] 10 15] ]
1220            }
1221            frac-set {
1222                if ![validreal value 6 4] {return 0}
1223                setexp ${key} $value 11 6
1224            }
1225            x-get {
1226                return [string trim [string range [readexp ${key}] 16 23] ]
1227            }
1228            x-set {
1229                if ![validreal value 8 5] {return 0}
1230                setexp ${key} $value 17 8
1231            }
1232            y-get {
1233                return [string trim [string range [readexp ${key}] 24 31] ]
1234            }
1235            y-set {
1236                if ![validreal value 8 5] {return 0}
1237                setexp ${key} $value 25 8
1238            }
1239            z-get {
1240                return [string trim [string range [readexp ${key}] 32 39] ]
1241            }
1242            z-set {
1243                if ![validreal value 8 5] {return 0}
1244                setexp ${key} $value 33 8
1245            }
1246            Uiso-get {
1247                return [string trim [string range [readexp ${key}] 40 45] ]
1248            }
1249            Uiso-set {
1250                if ![validreal value 6 4] {return 0}
1251                setexp ${key} $value 41 6
1252            }
1253            label-get {
1254                return [string trim [string range [readexp ${key}] 46 50] ]
1255            }
1256            label-set {
1257                setexp ${key} $value 47 5
1258            }
1259            residue-get {
1260                return [string range [readexp ${key}] 51 53]
1261            }
1262            residue-set {
1263                setexp ${key} $value 52 3
1264            }
1265            group-get {
1266                return [string range [readexp ${key}] 54 55]
1267            }
1268            group-set {
1269                setexp ${key} $value 55 2
1270            }
1271            resnum-get {
1272                return [string trim [string range [readexp ${key}] 56 59] ]
1273            }
1274            resnum-set {
1275                if ![validint value 4] {return 0}
1276                setexp "${key} EPHAS" $value 57 4
1277            }
1278            fref-get {
1279                if {[string toupper [string range [readexp $key] 60 60]] == "F"} {
1280                    return 1
1281                }
1282                return 0
1283            }
1284            fref-set {
1285                if $value {
1286                    setexp $key "F" 61 1
1287                } else {
1288                    setexp $key " " 61 1
1289                }           
1290            }
1291            xref-get {
1292                if {[string toupper [string range [readexp $key] 61 61]] == "X"} {
1293                    return 1
1294                }
1295                return 0
1296            }
1297            xref-set {
1298                if $value {
1299                    setexp $key "X" 62 1
1300                } else {
1301                    setexp ${key}B " " 62 1
1302                }           
1303            }
1304            uref-get {
1305                if {[string toupper [string range [readexp $key] 62 62]] == "U"} {
1306                    return 1
1307                }
1308                return 0
1309            }
1310            uref-set {
1311                if $value {
1312                    setexp $key "U" 63 1
1313                } else {
1314                    setexp $key " " 63 1
1315                }           
1316            }
1317
1318            fdamp-get {
1319                set val [string range [readexp ${key}] 63 63]
1320                if {$val == " "} {return 0}
1321                return $val
1322            }
1323            fdamp-set {
1324                setexp ${key} $value 64 1
1325            }
1326            xdamp-get {
1327                set val [string range [readexp ${key}] 64 64]
1328                if {$val == " "} {return 0}
1329                return $val
1330            }
1331            xdamp-set {
1332                setexp ${key} $value 65 1
1333            }
1334
1335            udamp-get {
1336                set val [string range [readexp ${key}] 65 65]
1337                if {$val == " "} {return 0}
1338                return $val
1339            }
1340            udamp-set {
1341                setexp ${key} $value 66 1
1342            }
1343            default {
1344                set msg "Unsupported mmatominfo access: parm=$parm action=$action"
1345                tk_dialog .badexp "Error in readexp" $msg error 0 Exit
1346            }
1347        }
1348    }
1349    return 1
1350}
1351
1352
1353
1354# get histogram information: histinfo histlist parm action value
1355# histlist is a list of histogram numbers
1356# parm:
1357#     title
1358#     file  -- file name of raw data for histogram (*)
1359#     scale (*)
1360#     sref/sdamp -- refinement flag/damping value for the scale factor (*)
1361#     lam1, lam2 (*)
1362#     ttref refinement flag for the 2theta (ED Xray) (*)
1363#     wref refinement flag for the wavelength (*)
1364#     ratref refinement flag for the wavelength ratio (*)
1365#     difc, difa -- TOF calibration constants (*)
1366#     dcref,daref -- refinement flag for difc, difa (*)
1367#     zero (*)
1368#     zref refinement flag for the zero correction (*)
1369#     ipola (*)
1370#     pola (*)
1371#     pref refinement flag for the polarization (*)
1372#     kratio (*)
1373#     ddamp -- damping value for the diffractometer constants (*)
1374#     backtype -- background function number *
1375#     backterms -- number of background terms *
1376#     bref/bdamp -- refinement flag/damping value for the background (*)
1377#     bterm$n -- background term #n (*)
1378#     bank -- Bank number
1379#     tofangle -- detector angle (TOF only)
1380#     foextract  -- Fobs extraction flag (*)
1381#     LBdamp  -- LeBail damping value (*)
1382#     tmin/tmax -- minimum & maximum usable 2theta/TOF/energy
1383#     excl -- excluded regions (*)
1384#     dmin -- minimum d-space for reflection generation (*)
1385#     use  -- use flag; 1 = use; 0 = do not use (*)
1386#     dstart -- dummy histogram starting tmin/emin/2theta (*)
1387#     dstep -- dummy histogram step size tmin/emin/2theta (*)
1388#     dpoints -- dummy histogram number of points (*)
1389#     dtype   -- dummy histogram type (CONST or SLOG)
1390#     abscor1 -- 1st absorption correction (*)
1391#     abscor2 -- 2nd absorption correction (*)
1392#     abstype -- absorption correction type (*)
1393#     absdamp -- damping for absorption refinement (*)
1394#     absref -- refinement damping for absorption refinement (*)
1395#     proftbl -- returns number of profile table terms
1396#     anomff -- returns a list of elements, f' and f"
1397#   parameters transferred from the instrument parameter file:
1398#     ITYP    -- returns the contents of the ITYP record
1399proc histinfo {histlist parm "action get" "value {}"} {
1400    global expgui
1401    foreach hist $histlist {
1402        if {$hist < 10} {
1403            set key "HST  $hist"
1404        } else {
1405            set key "HST $hist"
1406        }
1407        switch -glob ${parm}-$action {
1408            foextract-get {
1409                set line [readexp "${key} EPHAS"]
1410                # add a EPHAS if not exists
1411                if {$line == {}} {
1412                    makeexprec "${key} EPHAS"
1413                    # expedt defaults this to "F", but I think "T" is better
1414                    setexp "${key} EPHAS" "Y" 50 1
1415                    if $expgui(debug) {puts "Warning: creating a ${key} EPHAS record"}
1416                }
1417                if {[string toupper [string range $line 49 49]] == "T"} {
1418                    return 1
1419                }
1420                # the flag has changed to "Y/N" in the latest versions
1421                # of GSAS
1422                if {[string toupper [string range $line 49 49]] == "Y"} {
1423                    return 1
1424                }
1425                return 0
1426            }
1427            foextract-set {
1428                # the flag has changed to "Y/N" in the latest versions
1429                # of GSAS
1430                if $value {
1431                    setexp "${key} EPHAS" "Y" 50 1
1432                } else {
1433                    setexp "${key} EPHAS" "N" 50 1
1434                }
1435            }
1436            LBdamp-get {
1437                set v [string trim [string range [readexp "${key} EPHAS"] 54 54]]
1438                if {$v == ""} {return 0}
1439                return $v
1440            }
1441            LBdamp-set {
1442                if ![validint value 5] {return 0}
1443                setexp "${key} EPHAS" $value 51 5
1444            }
1445            title-get {
1446                return [string trim [readexp "${key}  HNAM"] ]
1447            }
1448            scale-get {
1449                return [string trim [string range [readexp ${key}HSCALE] 0 14]]
1450            }
1451            scale-set {
1452                if ![validreal value 15 6] {return 0}
1453                setexp ${key}HSCALE $value 1 15
1454            }
1455            sref-get {
1456                if {[string toupper [string range [readexp ${key}HSCALE] 19 19]] == "Y"} {
1457                    return 1
1458                }
1459                return 0
1460            }
1461            sref-set {
1462                if $value {
1463                    setexp ${key}HSCALE "Y" 20 1
1464                } else {
1465                    setexp ${key}HSCALE "N" 20 1
1466                }           
1467            }
1468            sdamp-get {
1469                set val [string range [readexp ${key}HSCALE] 24 24]
1470                if {$val == " "} {return 0}
1471                return $val
1472            }
1473            sdamp-set {
1474                setexp ${key}HSCALE $value 25 1
1475            }
1476
1477            difc-get -
1478            lam1-get {
1479                return [string trim [string range [readexp "${key} ICONS"] 0 9]]
1480            }
1481            difc-set -
1482            lam1-set {
1483                if ![validreal value 10 7] {return 0}
1484                setexp "${key} ICONS" $value 1 10
1485                # set the powpref warning (1 = suggested)
1486                catch {
1487                    global expgui
1488                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
1489                    set msg "Diffractometer constants" 
1490                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1491                        append expgui(needpowpref_why) "\t$msg were changed\n"
1492                    }
1493                }
1494            }
1495            difa-get -
1496            lam2-get {
1497                return [string trim [string range [readexp "${key} ICONS"] 10 19]]
1498            }
1499            difa-set -
1500            lam2-set {
1501                if ![validreal value 10 7] {return 0}
1502                setexp "${key} ICONS" $value 11 10
1503                # set the powpref warning (1 = suggested)
1504                catch {
1505                    global expgui
1506                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
1507                    set msg "Diffractometer constants" 
1508                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1509                        append expgui(needpowpref_why) "\t$msg were changed\n"
1510                    }
1511                }
1512            }
1513            zero-get {
1514                return [string trim [string range [readexp "${key} ICONS"] 20 29]]
1515            }
1516            zero-set {
1517                if ![validreal value 10 5] {return 0}
1518                setexp "${key} ICONS" $value 21 10
1519                # set the powpref warning (1 = suggested)
1520                catch {
1521                    global expgui
1522                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
1523                    set msg "Diffractometer constants" 
1524                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1525                        append expgui(needpowpref_why) "\t$msg were changed\n"
1526                    }
1527                }
1528            }
1529            ipola-get {
1530                return [string trim [string range [readexp "${key} ICONS"] 54 54]]
1531            }
1532            ipola-set {
1533                if ![validint value 1] {return 0}
1534                setexp "${key} ICONS" $value 55 1
1535            }
1536            pola-get {
1537                return [string trim [string range [readexp "${key} ICONS"] 40 49]]
1538            }
1539            pola-set {
1540                if ![validreal value 10 5] {return 0}
1541                setexp "${key} ICONS" $value 41 10
1542            }
1543            kratio-get {
1544                set val [string trim [string range [readexp "${key} ICONS"] 55 64]]
1545                if {$val == ""} {set val 0}
1546                # N.B. this code is used w/CW, where Kratio may not be 0.0
1547                set lam2 [string trim [string range [readexp "${key} ICONS"] 10 19]]
1548                if {$lam2 == ""} {set lam2 0}
1549                # Change kratio & flag the change (this is rather kludged)
1550                if {$val == 0 && $lam2 != 0} {
1551                    set val 0.5
1552                    validreal val 10 5
1553                    setexp "${key} ICONS" $val 56 10
1554                    catch {incr ::expgui(changed)}
1555                }
1556                return $val
1557            }
1558            kratio-set {
1559                if ![validreal value 10 5] {return 0}
1560                setexp "${key} ICONS" $value 56 10
1561            }
1562
1563            wref-get {
1564            #------------------------------------------------------
1565            # col 33: refine flag for lambda, difc, ratio and theta
1566            #------------------------------------------------------
1567                if {[string toupper [string range \
1568                        [readexp "${key} ICONS"] 32 32]] == "L"} {
1569                    return 1
1570                }
1571                return 0
1572            }
1573            wref-set {
1574                if $value {
1575                    setexp "${key} ICONS" "L" 33 1
1576                } else {
1577                    setexp "${key} ICONS" " " 33 1
1578                }           
1579            }
1580            ratref-get {
1581                if {[string toupper [string range \
1582                        [readexp "${key} ICONS"] 32 32]] == "R"} {
1583                    return 1
1584                }
1585                return 0
1586            }
1587            ratref-set {
1588                if $value {
1589                    setexp "${key} ICONS" "R" 33 1
1590                } else {
1591                    setexp "${key} ICONS" " " 33 1
1592                }           
1593            }
1594            dcref-get {
1595                if {[string toupper [string range \
1596                        [readexp "${key} ICONS"] 32 32]] == "C"} {
1597                    return 1
1598                }
1599                return 0
1600            }
1601            dcref-set {
1602                if $value {
1603                    setexp "${key} ICONS" "C" 33 1
1604                } else {
1605                    setexp "${key} ICONS" " " 33 1
1606                }           
1607            }
1608            ttref-get {
1609                if {[string toupper [string range \
1610                        [readexp "${key} ICONS"] 32 32]] == "T"} {
1611                    return 1
1612                }
1613                return 0
1614            }
1615            ttref-set {
1616                if $value {
1617                    setexp "${key} ICONS" "T" 33 1
1618                } else {
1619                    setexp "${key} ICONS" " " 33 1
1620                }           
1621            }
1622
1623
1624            pref-get {
1625            #------------------------------------------------------
1626            # col 34: refine flag for POLA & DIFA
1627            #------------------------------------------------------
1628                if {[string toupper [string range \
1629                        [readexp "${key} ICONS"] 33 33]] == "P"} {
1630                    return 1
1631                }
1632                return 0
1633            }
1634            pref-set {
1635                if $value {
1636                    setexp "${key} ICONS" "P" 34 1
1637                } else {
1638                    setexp "${key} ICONS" " " 34 1
1639                }           
1640            }
1641            daref-get {
1642                if {[string toupper [string range \
1643                        [readexp "${key} ICONS"] 33 33]] == "A"} {
1644                    return 1
1645                }
1646                return 0
1647            }
1648            daref-set {
1649                if $value {
1650                    setexp "${key} ICONS" "A" 34 1
1651                } else {
1652                    setexp "${key} ICONS" " " 34 1
1653                }           
1654            }
1655
1656            zref-get {
1657            #------------------------------------------------------
1658            # col 34: refine flag for zero correction
1659            #------------------------------------------------------
1660                if {[string toupper [string range [readexp "${key} ICONS"] 34 34]] == "Z"} {
1661                    return 1
1662                }
1663                return 0
1664            }
1665            zref-set {
1666                if $value {
1667                    setexp "${key} ICONS" "Z" 35 1
1668                } else {
1669                    setexp "${key} ICONS" " " 35 1
1670                }           
1671            }
1672
1673            ddamp-get {
1674                set val [string range [readexp "${key} ICONS"] 39 39]
1675                if {$val == " "} {return 0}
1676                return $val
1677            }
1678            ddamp-set {
1679                setexp "${key} ICONS" $value 40 1
1680            }
1681
1682            backtype-get {
1683                set val [string trim [string range [readexp "${key}BAKGD "] 0 4]]
1684                if {$val == " "} {return 0}
1685                return $val
1686            }
1687            backtype-set {
1688                if ![validint value 5] {return 0}
1689                setexp "${key}BAKGD " $value 1 5
1690            }
1691            backterms-get {
1692                set val [string trim [string range [readexp "${key}BAKGD "] 5 9]]
1693                if {$val == " "} {return 0}
1694                return $val
1695            }
1696            backterms-set {
1697                # this takes a bit of work -- if terms are added, add lines as needed to the .EXP
1698                set oldval [string trim [string range [readexp "${key}BAKGD "] 5 9]]
1699                if ![validint value 5] {return 0}
1700                if {$oldval < $value} {
1701                    set line1  [expr {2 + ($oldval - 1) / 4}]
1702                    set line2  [expr {1 + ($value - 1) / 4}]
1703                    for {set i $line1} {$i <= $line2} {incr i} {
1704                        # create a blank entry if needed
1705                        makeexprec ${key}BAKGD$i
1706                    }
1707                    incr oldval
1708                    for {set num $oldval} {$num <= $value} {incr num} {
1709                        set f1 [expr {15*(($num - 1) % 4)}]
1710                        set f2 [expr {15*(1 + ($num - 1) % 4)-1}]
1711                        set line  [expr {1 + ($num - 1) / 4}]
1712                        if {[string trim [string range [readexp ${key}BAKGD$line] $f1 $f2]] == ""} {
1713                            set f1 [expr {15*(($num - 1) % 4)+1}]
1714                            setexp ${key}BAKGD$line 0.0 $f1 15                 
1715                        }
1716                    }
1717                }
1718                setexp "${key}BAKGD " $value 6 5
1719
1720            }
1721            bref-get {
1722                if {[string toupper [string range [readexp "${key}BAKGD"] 14 14]] == "Y"} {
1723                    return 1
1724                }
1725                return 0
1726            }
1727            bref-set {
1728                if $value {
1729                    setexp "${key}BAKGD "  "Y" 15 1
1730                } else {
1731                    setexp "${key}BAKGD "  "N" 15 1
1732                }
1733            }
1734            bdamp-get {
1735                set val [string range [readexp "${key}BAKGD "] 19 19]
1736                if {$val == " "} {return 0}
1737                return $val
1738            }
1739            bdamp-set {
1740                setexp "${key}BAKGD " $value 20 1
1741            }
1742            bterm*-get {
1743                regsub bterm $parm {} num
1744                set f1 [expr {15*(($num - 1) % 4)}]
1745                set f2 [expr {15*(1 + ($num - 1) % 4)-1}]
1746                set line  [expr {1 + ($num - 1) / 4}]
1747                return [string trim [string range [readexp ${key}BAKGD$line] $f1 $f2] ]
1748            }
1749            bterm*-set {
1750                regsub bterm $parm {} num
1751                if ![validreal value 15 6] {return 0}
1752                set f1 [expr {15*(($num - 1) % 4)+1}]
1753                set line  [expr {1 + ($num - 1) / 4}]
1754                setexp ${key}BAKGD$line $value $f1 15
1755            }
1756            bank-get {
1757                return [string trim [string range [readexp "${key} BANK"] 0 4]]
1758            }
1759            tofangle-get {
1760                return [string trim [string range [readexp "${key}BNKPAR"] 10 19]]
1761            }
1762            tmin-get {
1763                return [string trim [string range [readexp "${key} TRNGE"] 0 9]]
1764            }
1765            tmax-get {
1766                return [string trim [string range [readexp "${key} TRNGE"] 10 19]]
1767            }
1768            excl-get {
1769                set n [string trim [string range [readexp "${key} NEXC"] 0 4]]
1770                set exlist {}
1771                for {set i 1} {$i <= $n} {incr i} {
1772                    set line [readexp [format "${key}EXC%3d" $i]]
1773                    lappend exlist [list \
1774                            [string trim [string range $line  0  9]] \
1775                            [string trim [string range $line 10 19]]]
1776                }
1777                return $exlist
1778            }
1779            excl-set {
1780                set n [llength $value]
1781                if ![validint n 5] {return 0}
1782                setexp "${key} NEXC" $n 1 5
1783                set i 0
1784                foreach p $value {
1785                    incr i
1786                    foreach {r1 r2} $p {}
1787                    validreal r1 10 3
1788                    validreal r2 10 3
1789                    set k [format "${key}EXC%3d" $i]
1790                    if {![existsexp $k]} {
1791                        makeexprec $k
1792                    }
1793                    setexp $k ${r1}${r2} 1 20
1794                }
1795                # set the powpref warning (2 = required)
1796                catch {
1797                    global expgui
1798                    set expgui(needpowpref) 2
1799                    set msg "Excluded regions" 
1800                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1801                        append expgui(needpowpref_why) "\t$msg were changed\n"
1802                    }
1803                }
1804            }
1805            file-get {
1806                return [string trim [readexp "${key}  HFIL"] ]
1807            }
1808            file-set {
1809                setexp "${key}  HFIL" $value 3 65
1810            }
1811            bank-get {
1812                return [string trim [string range [readexp "${key} BANK"] 0 4]]
1813            }
1814            dmin-get {
1815                return [string trim [string range [readexp "${key} NREF"] 5 14]]
1816            }
1817            dmin-set {
1818                if ![validreal value 10 4] {return 0}
1819                setexp "${key} NREF" $value 6 10
1820                # set the powpref warning (2 = required)
1821                catch {
1822                    global expgui
1823                    set expgui(needpowpref) 2
1824                    set msg "Dmin (reflection range)" 
1825                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1826                        append expgui(needpowpref_why) "\t$msg was changed\n"
1827                    }
1828                }
1829            }
1830            use-get {
1831                set k [expr {($hist+11)/12}]
1832                set line [readexp " EXPR  HTYP$k"]
1833                set j [expr {((($hist-1) % 12)+1)*5}]
1834                if {[string range $line $j $j] == "*"} {return 0}
1835                return 1
1836            }
1837            use-set {
1838                set k [expr {($hist+11)/12}]
1839                set line [readexp " EXPR  HTYP$k"]
1840                set j [expr {((($hist-1) % 12)+1)*5+1}]
1841                if {$value} {
1842                    setexp " EXPR  HTYP$k" " " $j 1
1843                } else {
1844                    setexp " EXPR  HTYP$k" "*" $j 1
1845                }
1846                # set the powpref warning (2 = required)
1847                catch {
1848                    global expgui
1849                    set expgui(needpowpref) 2
1850                    set msg "Histogram use flags" 
1851                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1852                        append expgui(needpowpref_why) "\t$msg were changed\n"
1853                    }
1854                }
1855            }
1856            dstart-get {
1857                return [string trim [string range [readexp "${key} DUMMY"] 20 29]]
1858            }
1859            dstart-set {
1860                if ![validreal value 10 3] {return 0}
1861                setexp "${key} DUMMY" $value 21 10
1862                # set the powpref warning (1 = suggested)
1863                catch {
1864                    global expgui
1865                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
1866                    set msg "Dummy histogram parameters" 
1867                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1868                        append expgui(needpowpref_why) "\t$msg were changed\n"
1869                    }
1870                }
1871            }
1872            dstep-get {
1873                return [string trim [string range [readexp "${key} DUMMY"] 30 39]]
1874            }
1875            dstep-set {
1876                if ![validreal value 10 3] {return 0}
1877                setexp "${key} DUMMY" $value 31 10
1878                catch {
1879                    global expgui
1880                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
1881                    set msg "Dummy histogram parameters" 
1882                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1883                        append expgui(needpowpref_why) "\t$msg were changed\n"
1884                    }
1885                }
1886            }
1887            dpoints-get {
1888                return [string trim [string range [readexp "${key} DUMMY"] 0 9]]
1889            }
1890            dpoints-set {
1891                if ![validint value 10] {return 0}
1892                setexp "${key} DUMMY" $value 1 10
1893                catch {
1894                    global expgui
1895                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
1896                    set msg "Dummy histogram parameters" 
1897                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
1898                        append expgui(needpowpref_why) "\t$msg were changed\n"
1899                    }
1900                }
1901            }
1902            dtype-get {
1903                return [string trim [string range [readexp "${key} DUMMY"] 10 19]]
1904            }
1905            abscor1-get {
1906                return [string trim [string range [readexp "${key}ABSCOR"] 0 14]]
1907            }
1908            abscor1-set {
1909                if ![validreal value 15 7] {return 0}
1910                setexp "${key}ABSCOR" $value 1 15
1911            }
1912            abscor2-get {
1913                return [string trim [string range [readexp "${key}ABSCOR"] 15 29]]
1914            }
1915            abscor2-set {
1916                # this must have a decimal as the 5th character, so that we end up with a
1917                # decimal point in column 20.
1918                set tmp $value
1919                if ![validreal tmp 12 7] {return 0}
1920                set pos [string first "." $tmp]
1921                while {$pos < 4} {
1922                    set tmp " $tmp"
1923                    set pos [string first "." $tmp]
1924                }
1925                if {$pos == 4} {
1926                    setexp "${key}ABSCOR" $tmp 16 15
1927                    return
1928                }
1929                catch {
1930                    set tmp [format "%12.6E" $value]
1931                    set pos [string first "." $tmp]
1932                    while {$pos < 4} {
1933                        set tmp " $tmp"
1934                        set pos [string first "." $tmp]
1935                    }
1936                    if {$pos == 4} {
1937                        setexp "${key}ABSCOR" $tmp 16 15
1938                        return
1939                    }
1940                }
1941                return 0
1942            }
1943            abstype-get {
1944                set val [string trim [string range [readexp "${key}ABSCOR"] 40 44]]
1945                if {$val == ""} {set val 0}
1946                return $val
1947            }
1948            abstype-set {
1949                if ![validint value 5] {return 0}
1950                setexp "${key}ABSCOR" $value 41 5
1951            }
1952            absdamp-get {
1953                set val [string range [readexp "${key}ABSCOR"] 39 39]
1954                if {$val == " "} {return 0}
1955                return $val
1956            }
1957            absdamp-set {
1958                if ![validint value 5] {return 0}
1959                setexp "${key}ABSCOR" $value 36 5
1960            }
1961            absref-get {
1962                if {[string toupper \
1963                        [string range [readexp "${key}ABSCOR"] 34 34]] == "Y"} {
1964                    return 1
1965                }
1966                return 0
1967            }
1968            absref-set {
1969                if $value {
1970                    setexp "${key}ABSCOR" "    Y" 31 5
1971                } else {
1972                    setexp "${key}ABSCOR" "    N" 31 5
1973                }
1974            }
1975            ITYP-get {
1976                return [string trim [readexp "${key}I ITYP"]]
1977            }
1978            proftbl-get {
1979                set line [readexp "${key}PAB3"]
1980                if {$line == ""} {return 0}
1981                set val [string trim [string range $line 0 4]]
1982                if {$val == ""} {return 0}
1983                return $val
1984            }
1985            anomff-get {
1986                set l {}
1987                foreach i {1 2 3 4 5 6 7 8 9} {
1988                    if {![existsexp "${key}FFANS$i"]} continue
1989                    set line [readexp "${key}FFANS$i"]
1990                    set elem [string trim [string range $line 2 9]]
1991                    set fp [string trim [string range $line 10 19]]
1992                    set fpp [string trim [string range $line 20 29]]
1993                    lappend l [list $elem $fp $fpp]
1994                }
1995                return $l
1996            }
1997            anomff-set {
1998                # match up input against elements in list.
1999                # change elements included, return any elements that are
2000                # not found.
2001                set errorlist {}
2002                foreach triplet $value {
2003                    foreach {e fp fpp} $triplet {}               
2004                    foreach i {1 2 3 4 5 6 7 8 9} {
2005                        if {![existsexp "${key}FFANS$i"]} continue
2006                        # note that the element name is not used or validated
2007                        set elem [string trim [string range \
2008                                                   [readexp "${key}FFANS$i"] 2 9]]
2009                        if {[string match -nocase $e $elem]} { 
2010                            if ![validreal fp 10 3] {return 0}
2011                            setexp "${key}FFANS$i" $fp 11 10
2012                            if ![validreal fpp 10 3] {return 0}
2013                            setexp "${key}FFANS$i" $fpp 21 10
2014                            set e {}
2015                            break
2016                        }
2017                    }
2018                    if {$e != ""} {
2019                        # oops, no match
2020                        lappend errorlist $e
2021                    }
2022                }
2023                if {$errorlist != ""} {return [list 0 $errorlist]}
2024            }
2025            default {
2026                set msg "Unsupported histinfo access: parm=$parm action=$action"
2027                tk_dialog .badexp "Error in readexp" $msg error 0 Exit
2028            }
2029        }
2030    }
2031    return 1
2032}
2033
2034# read the information that differs by both histogram and phase (profile & phase fraction)
2035# use: hapinfo hist phase parm action value
2036
2037#     frac -- phase fraction (*)
2038#     frref/frdamp -- refinement flag/damping value for the phase fraction (*)
2039#     proftype -- profile function number (*)
2040#     profterms -- number of profile terms (*)
2041#     pdamp -- damping value for the profile (*)
2042#     pcut -- cutoff value for the profile (*)
2043#     pterm$n -- profile term #n (*)
2044#     pref$n -- refinement flag value for profile term #n (*)
2045#     extmeth -- Fobs extraction method (*)
2046#     POnaxis -- number of defined M-D preferred axes
2047proc hapinfo {histlist phaselist parm "action get" "value {}"} {
2048    foreach phase $phaselist hist $histlist {
2049        if {$phase == ""} {set phase [lindex $phaselist end]}
2050        if {$hist == ""} {set hist [lindex $histlist end]}
2051        if {$hist < 10} {
2052            set hist " $hist"
2053        }
2054        set key "HAP${phase}${hist}"
2055        switch -glob ${parm}-$action {
2056            extmeth-get {
2057                set i1 [expr {($phase - 1)*5}]
2058                set i2 [expr {$i1 + 4}]
2059                return [string trim [string range [readexp "HST $hist EPHAS"] $i1 $i2]]
2060            }
2061            extmeth-set {
2062                set i1 [expr {($phase - 1)*5 + 1}]
2063                if ![validint value 5] {return 0}
2064                setexp "HST $hist EPHAS" $value $i1 5
2065            }
2066            frac-get {
2067                return [string trim [string range [readexp ${key}PHSFR] 0 14]]
2068            }
2069            frac-set {
2070                if ![validreal value 15 6] {return 0}
2071                setexp ${key}PHSFR $value 1 15
2072            }
2073            frref-get {
2074                if {[string toupper [string range [readexp ${key}PHSFR] 19 19]] == "Y"} {
2075                    return 1
2076                }
2077                return 0
2078            }
2079            frref-set {
2080                if $value {
2081                    setexp ${key}PHSFR "Y" 20 1
2082                } else {
2083                    setexp ${key}PHSFR "N" 20 1
2084                }           
2085            }
2086            frdamp-get {
2087                set val [string range [readexp ${key}PHSFR] 24 24]
2088                if {$val == " "} {return 0}
2089                return $val
2090            }
2091            frdamp-set {
2092                setexp ${key}PHSFR $value 25 1
2093            }
2094            proftype-get {
2095                set val [string range [readexp "${key}PRCF "] 0 4]
2096                if {$val == " "} {return 0}
2097                return $val
2098            }
2099            proftype-set {
2100                if ![validint value 5] {return 0}
2101                setexp "${key}PRCF " $value 1 5
2102                # set the powpref warning (1 = suggested)
2103                catch {
2104                    global expgui
2105                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
2106                    set msg "Profile parameters" 
2107                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
2108                        append expgui(needpowpref_why) "\t$msg were changed\n"
2109                    }
2110                }
2111            }
2112            profterms-get {
2113                set val [string range [readexp "${key}PRCF "] 5 9]
2114                if {$val == " "} {return 0}
2115                return $val
2116            }
2117            profterms-set {
2118                if ![validint value 5] {return 0}
2119                setexp "${key}PRCF " $value 6 5
2120                # now check that all needed entries exist
2121                set lines [expr {1 + ($value - 1) / 4}]
2122                for {set i 1} {$i <= $lines} {incr i} {
2123                    makeexprec "${key}PRCF $i"
2124                }
2125                # set the powpref warning (1 = suggested)
2126                catch {
2127                    global expgui
2128                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
2129                    set msg "Profile parameters" 
2130                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
2131                        append expgui(needpowpref_why) "\t$msg were changed\n"
2132                    }
2133                }
2134            }
2135            pcut-get {
2136                return [string trim [string range [readexp "${key}PRCF "] 10 19]]
2137            }
2138            pcut-set {
2139                if ![validreal value 10 5] {return 0}
2140                setexp "${key}PRCF " $value 11 10
2141                # set the powpref warning (1 = suggested)
2142                catch {
2143                    global expgui
2144                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
2145                    set msg "Profile parameters" 
2146                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
2147                        append expgui(needpowpref_why) "\t$msg were changed\n"
2148                    }
2149                }
2150            }
2151            pdamp-get {
2152                set val [string range [readexp "${key}PRCF "] 24 24]
2153                if {$val == " "} {return 0}
2154                return $val
2155            }
2156            pdamp-set {
2157                setexp "${key}PRCF   " $value 25 1
2158            }
2159            pterm*-get {
2160                regsub pterm $parm {} num
2161                set f1 [expr {15*(($num - 1) % 4)}]
2162                set f2 [expr {15*(1 + ($num - 1) % 4)-1}]
2163                set line  [expr {1 + ($num - 1) / 4}]
2164                return [string trim [string range [readexp "${key}PRCF $line"] $f1 $f2] ]
2165            }
2166            pterm*-set {
2167                if ![validreal value 15 6] {return 0}
2168                regsub pterm $parm {} num
2169                set f1 [expr {1+ 15*(($num - 1) % 4)}]
2170                set line  [expr {1 + ($num - 1) / 4}]
2171                setexp "${key}PRCF $line" $value $f1 15
2172                # set the powpref warning (1 = suggested)
2173                catch {
2174                    global expgui
2175                    if {$expgui(needpowpref) == 0} {set expgui(needpowpref) 1}
2176                    set msg "Profile parameters" 
2177                    if {[string first $msg $expgui(needpowpref_why)] == -1} {
2178                        append expgui(needpowpref_why) "\t$msg were changed\n"
2179                    }
2180                }
2181            }
2182            pref*-get {
2183                regsub pref $parm {} num
2184                set f [expr {24+$num}]
2185                if {[string toupper [string range [readexp "${key}PRCF  "] $f $f]] == "Y"} {
2186                    return 1
2187                }
2188                return 0
2189            }
2190            pref*-set {
2191                regsub pref $parm {} num
2192                set f [expr {25+$num}]
2193                if $value {
2194                    setexp ${key}PRCF "Y" $f 1
2195                } else {
2196                    setexp ${key}PRCF "N" $f 1
2197                }           
2198            }
2199            POnaxis-get {
2200                set val [string trim \
2201                        [string range [readexp "${key}NAXIS"] 0 4]]
2202                if {$val == ""} {return 0}
2203                return $val
2204            }
2205            POnaxis-set {
2206                if ![validint value 5] {return 0}
2207                # there should be a NAXIS record, but if not make one
2208                if {![existsexp "${key}NAXIS"]} {
2209                    makeexprec "${key}NAXIS"
2210                }
2211                setexp "${key}NAXIS  " $value 1 5
2212            }
2213            default {
2214                set msg "Unsupported hapinfo access: parm=$parm action=$action"
2215                tk_dialog .badexp "Error in readexp" $msg error 0 Exit
2216            }
2217        }
2218    }
2219    return 1
2220}
2221
2222#  read fixed constraints
2223
2224proc atom_constraint_read {phase} {
2225    set fix_list ""
2226    foreach k {1 2 3 4 5 6 7 8 9} {
2227        set key [format "LEQV HOLD%1d%2d" $phase $k]
2228        set line [readexp $key]
2229        foreach j {2 10 18 26 34 42 50 58} {
2230            set fix_param [string range $line $j [expr $j+7]]
2231            if {[string trim $fix_param] == ""} {return $fix_list}
2232            lappend fix_list $fix_param
2233        }
2234    }
2235}
2236
2237# load all atom constraints into global array fix_param
2238proc atom_constraint_load { } {
2239    catch {unset ::fix_param}
2240    foreach i $::expmap(phaselist) {
2241        set temp [atom_constraint_read $i]
2242        foreach j $temp {
2243            set atomnum [string trim [string range $j 2 3]]
2244            set param [string trim [string range $j 4 6]]
2245            set ::fix_param($i,$atomnum,$param) 1   
2246        }
2247    }
2248}
2249
2250proc atom_constraint_write {phase fix_list} {
2251    foreach key [array names ::exparray "LEQV HOLD$phase*"] {
2252        delexp $key
2253    }
2254    set k 0
2255    set j 1
2256    set line ""
2257    foreach fix $fix_list {
2258        incr k
2259        append line $fix
2260        if {$k == 8} {
2261            set key [format "LEQV HOLD%1d%2d" $phase $j]
2262            makeexprec $key
2263            setexp $key $line 3 [expr ($k * 8) + 2]
2264            set k 0
2265            incr j
2266            set line ""
2267        }
2268    }
2269    if {$line != ""} {
2270        set key [format "LEQV HOLD%1d%2d" $phase $j]
2271        makeexprec $key
2272        setexp $key $line 3 [expr ($k * 8) + 2]
2273    }   
2274}
2275
2276
2277#  get a logical constraint
2278#
2279#  type action
2280#  -----------
2281#  atom get  number        returns a list of constraints.
2282#   "   set  number value  replaces a list of constraints
2283#                          (value is a list of constraints)
2284#   "   add  number value  inserts a new list of constraints
2285#                          (number is ignored)
2286#   "   delete number      deletes a set of constraint entries
2287# Each item in the list of constraints is composed of 4 items:
2288#              phase, atom, variable, multiplier
2289# If variable=UISO atom can be ALL, otherwise atom is a number
2290# legal variable names: FRAC, X, Y, Z, UISO, U11, U22, U33, U12, U23, U13,
2291#                       MX, MY, MZ
2292#
2293#  type action
2294#  -----------
2295#  profileXX get number         returns a list of constraints for term XX=1-36
2296#                               use number=0 to get # of defined
2297#                                  constraints for term XX
2298#   "        set number value   replaces a list of constraints
2299#                               (value is a list of constraints)
2300#   "        add number value   inserts a new list of constraints
2301#                               (number is ignored)
2302#   "        delete number      deletes a set of constraint entries
2303# Each item in the list of constraints is composed of 3 items:
2304#              phase-list, histogram-list, multiplier
2305# Note that phase-list and/or histogram-list can be ALL
2306
2307proc constrinfo {type action number "value {}"} {
2308    global expmap
2309    if {[lindex $expmap(phasetype) 0] == 4} {
2310        set mm 1
2311    } else {
2312        set mm 0
2313    }
2314    switch -glob ${type}-$action {
2315        atom-get {
2316            # does this constraint exist?
2317            set key [format "LNCN%4d%4d" $number 1]
2318            if {![existsexp $key]} {return -1}
2319            set clist {}
2320            for {set i 1} {$i < 999} {incr i} {
2321                set key [format "LNCN%4d%4d" $number $i]
2322                if {![existsexp $key]} break
2323                set line [readexp $key]
2324                set j1 2
2325                set j2 17
2326                set seg [string range $line $j1 $j2]
2327                while {[string trim $seg] != ""} {
2328                    set p [string range $seg 0 0]
2329                    if {$p == 1 && $mm} {
2330                        set atom [string trim [string range $seg 1 4]]
2331                        set var [string trim [string range $seg 5 7]]
2332                        if {$atom == "ALL"} {
2333                            set var UIS
2334                        } else {
2335                            scan $atom %x atom
2336                        }
2337                        lappend clist [list $p $atom $var \
2338                                [string trim [string range $seg 8 end]]]
2339                    } else {
2340                        lappend clist [list $p \
2341                                [string trim [string range $seg 1 3]] \
2342                                [string trim [string range $seg 4 7]] \
2343                                [string trim [string range $seg 8 end]]]
2344                    }
2345                    incr j1 16
2346                    incr j2 16
2347                    set seg [string range $line $j1 $j2]
2348                }
2349            }
2350            return $clist
2351        }
2352        atom-set {
2353            # delete records for current constraint
2354            for {set i 1} {$i < 999} {incr i} {
2355                set key [format "LNCN%4d%4d" $number $i]
2356                if {![existsexp $key]} break
2357                delexp $key
2358            }
2359            set line {}
2360            set i 1
2361            foreach tuple $value {
2362                set p [lindex $tuple 0]
2363                if {$p == 1 && $mm && \
2364                        [string toupper [lindex $tuple 1]] == "ALL"} {
2365                    set seg [format %1dALL UIS%8.4f \
2366                            [lindex $tuple 0] \
2367                            [lindex $tuple 3]]
2368                } elseif {$p == 1 && $mm} {
2369                    set seg [eval format %1d%.4X%-3s%8.4f $tuple]
2370                } elseif {[string toupper [lindex $tuple 1]] == "ALL"} {
2371                    set seg [format %1dALL%-4s%8.4f \
2372                            [lindex $tuple 0] \
2373                            [lindex $tuple 2] \
2374                            [lindex $tuple 3]]
2375                } else {
2376                    set seg [eval format %1d%3d%-4s%8.4f $tuple]
2377                }
2378                append line $seg
2379                if {[string length $line] > 50} {
2380                    set key  [format "LNCN%4d%4d" $number $i]
2381                    makeexprec $key
2382                    setexp $key $line 3 68
2383                    set line {}
2384                    incr i
2385                }
2386            }
2387            if {$line != ""} {
2388                set key  [format "LNCN%4d%4d" $number $i]
2389                makeexprec $key
2390                setexp $key $line 3 68
2391            }
2392            return
2393        }
2394        atom-add {
2395            # loop over defined constraints
2396            for {set j 1} {$j < 9999} {incr j} {
2397                set key [format "LNCN%4d%4d" $j 1]
2398                if {![existsexp $key]} break
2399            }
2400            set number $j
2401            # save the constraint
2402            set line {}
2403            set i 1
2404            foreach tuple $value {
2405                set p [lindex $tuple 0]
2406                if {$p == 1 && $mm && \
2407                        [string toupper [lindex $tuple 1]] == "ALL"} {
2408                    set seg [format %1dALL UIS%8.4f \
2409                            [lindex $tuple 0] \
2410                            [lindex $tuple 3]]
2411                } elseif {$p == 1 && $mm} {
2412                    set seg [eval format %1d%.4X%-3s%8.4f $tuple]
2413                } elseif {[string toupper [lindex $tuple 1]] == "ALL"} {
2414                    set seg [format %1dALL%-4s%8.4f \
2415                            [lindex $tuple 0] \
2416                            [lindex $tuple 2] \
2417                            [lindex $tuple 3]]
2418                } else {
2419                    set seg [eval format %1d%3d%-4s%8.4f $tuple]
2420                }
2421                append line $seg
2422                if {[string length $line] > 50} {
2423                    set key  [format "LNCN%4d%4d" $number $i]
2424                    makeexprec $key
2425                    setexp $key $line 3 68
2426                    set line {}
2427                    incr i
2428                }
2429            }
2430            if {$line != ""} {
2431                set key  [format "LNCN%4d%4d" $number $i]
2432                makeexprec $key
2433                setexp $key $line 3 68
2434            }
2435            return
2436        }
2437        atom-delete {
2438            for {set j $number} {$j < 9999} {incr j} {
2439                # delete records for current constraint
2440                for {set i 1} {$i < 999} {incr i} {
2441                    set key [format "LNCN%4d%4d" $j $i]
2442                    if {![existsexp $key]} break
2443                    delexp $key
2444                }
2445                # now copy records, from the next entry, if any
2446                set j1 $j
2447                incr j1
2448                set key1 [format "LNCN%4d%4d" $j1 1]
2449                # if there is no record, there is nothing to copy -- done
2450                if {![existsexp $key1]} return
2451                for {set i 1} {$i < 999} {incr i} {
2452                    set key1 [format "LNCN%4d%4d" $j1 $i]
2453                    if {![existsexp $key1]} break
2454                    set key  [format "LNCN%4d%4d" $j  $i]
2455                    makeexprec $key
2456                    setexp $key [readexp $key1] 1 68
2457                }
2458            }
2459        }
2460        profile*-delete {
2461            regsub profile $type {} term
2462            if {$term < 10} {
2463                set term " $term"
2464            }
2465            set key "LEQV PF$term   "
2466            # return nothing if no term exists
2467            if {![existsexp $key]} {return 0}
2468
2469            # number of constraint terms
2470            set nterms [string trim [string range [readexp ${key}] 0 4] ]
2471            # don't delete a non-existing entry
2472            if {$number > $nterms} {return 0}
2473            set val [expr {$nterms - 1}]
2474            validint val 5
2475            setexp $key $val 1 5
2476            for {set i1 $number} {$i1 < $nterms} {incr i1} {
2477                set i2 [expr {1 + $i1}]
2478                # move the contents of constraint #i2 -> i1
2479                if {$i1 > 9} {
2480                    set k1 [expr {($i1+1)/10}]
2481                    set l1 $i1
2482                } else {
2483                    set k1 " "
2484                    set l1 " $i1"
2485                }
2486                set key1 "LEQV PF$term  $k1"
2487                # number of constraint lines for #i1
2488                set n1 [string trim [string range [readexp ${key1}] \
2489                        [expr {($i1%10)*5}] [expr {4+(($i1%10)*5)}]] ]
2490                if {$i2 > 9} {
2491                    set k2 [expr {($i2+1)/10}]
2492                    set l2 $i2
2493                } else {
2494                    set k2 " "
2495                    set l2 " $i2"
2496                }
2497                set key2 "LEQV PF$term  $k2"
2498                # number of constraint lines for #i2
2499                set n2 [string trim [string range [readexp ${key2}] \
2500                        [expr {($i2%10)*5}] [expr {4+(($i2%10)*5)}]] ]
2501                set val $n2
2502                validint val 5
2503                # move the # of terms
2504                setexp $key1 $val [expr {1+(($i1%10)*5)}] 5
2505                # move the terms
2506                for {set j 1} {$j <= $n2} {incr j 1} {
2507                    set key "LEQV PF${term}${l1}$j"
2508                    makeexprec $key
2509                    setexp $key [readexp "LEQV PF${term}${l2}$j"] 1 68
2510                }
2511                # delete any remaining lines
2512                for {set j [expr {$n2+1}]} {$j <= $n1} {incr j 1} {
2513                    delexp "LEQV PF${term}${l1}$j"
2514                }
2515            }
2516
2517            # clear the last term
2518            if {$nterms > 9} {
2519                set i [expr {($nterms+1)/10}]
2520            } else {
2521                set i " "
2522            }
2523            set key "LEQV PF$term  $i"
2524            set cb [expr {($nterms%10)*5}]
2525            set ce [expr {4+(($nterms%10)*5)}]
2526            set n2 [string trim [string range [readexp ${key}] $cb $ce] ]
2527            incr cb
2528            setexp $key "     " $cb 5
2529            # delete any remaining lines
2530            for {set j 1} {$j <= $n2} {incr j 1} {
2531                delexp "LEQV PF${term}${nterms}$j"
2532            }
2533        }
2534        profile*-set {
2535            regsub profile $type {} term
2536            if {$term < 10} {
2537                set term " $term"
2538            }
2539            set key "LEQV PF$term   "
2540            # get number of constraint terms
2541            set nterms [string trim [string range [readexp ${key}] 0 4] ]
2542            # don't change a non-existing entry
2543            if {$number > $nterms} {return 0}
2544            if {$number > 9} {
2545                set k1 [expr {($number+1)/10}]
2546                set l1 $number
2547            } else {
2548                set k1 " "
2549                set l1 " $number"
2550            }
2551            set key1 "LEQV PF$term  $k1"
2552            # old number of constraint lines
2553            set n1 [string trim [string range [readexp ${key1}] \
2554                    [expr {($number%10)*5}] [expr {4+(($number%10)*5)}]] ]
2555            # number of new constraints
2556            set j2 [llength $value]
2557            # number of new constraint lines
2558            set val [set n2 [expr {($j2 + 2)/3}]]
2559            # store the new # of lines
2560            validint val 5
2561            setexp $key1 $val [expr {1+(($number%10)*5)}] 5
2562
2563            # loop over the # of lines in the old or new, whichever is greater
2564            set v0 0
2565            for {set j 1} {$j <= [expr {($n1 > $n2) ? $n1 : $n2}]} {incr j 1} {
2566                set key "LEQV PF${term}${l1}$j"
2567                # were there more lines in the old?
2568                if {$j > $n2} {
2569                    # this line is not needed
2570                    if {$j % 3 == 1} {
2571                        delexp %key
2572                    }
2573                    continue
2574                }
2575                # are we adding new lines?
2576                if {$j > $n1} {
2577                    makeexprec $key
2578                }
2579                # add the three constraints to the line
2580                foreach s {3 23 43} \
2581                        item [lrange $value $v0 [expr {2+$v0}]] {
2582                    if {$item != ""} {
2583                        set val [format %-10s%9.3f \
2584                                [lindex $item 0],[lindex $item 1] \
2585                                [lindex $item 2]]
2586                        setexp $key $val $s 19
2587                    } else {
2588                        setexp $key " " $s 19
2589                    }
2590                }
2591                incr v0 3
2592            }
2593        }
2594        profile*-add {
2595            regsub profile $type {} term
2596            if {$term < 10} {
2597                set term " $term"
2598            }
2599            set key "LEQV PF$term   "
2600            if {![existsexp $key]} {makeexprec $key}
2601            set nterms [string trim [string range [readexp ${key}] 0 4] ]
2602            if {$nterms == ""} {
2603                set nterms 1
2604            } elseif {$nterms >= 99} {
2605                return 0
2606            } else {
2607                incr nterms
2608            }
2609            # store the new # of constraints
2610            set val $nterms
2611            validint val 5
2612            setexp $key $val 1 5
2613
2614            if {$nterms > 9} {
2615                set k1 [expr {($nterms+1)/10}]
2616                set l1 $nterms
2617            } else {
2618                set k1 " "
2619                set l1 " $nterms"
2620            }
2621            set key1 "LEQV PF$term  $k1"
2622
2623            # number of new constraints
2624            set j2 [llength $value]
2625            # number of new constraint lines
2626            set val [set n2 [expr {($j2 + 2)/3}]]
2627            # store the new # of lines
2628            validint val 5
2629            setexp $key1 $val [expr {1+(($nterms%10)*5)}] 5
2630
2631            # loop over the # of lines to be added
2632            set v0 0
2633            for {set j 1} {$j <= $n2} {incr j 1} {
2634                set key "LEQV PF${term}${l1}$j"
2635                makeexprec $key
2636                # add the three constraints to the line
2637                foreach s {3 23 43} \
2638                        item [lrange $value $v0 [expr {2+$v0}]] {
2639                    if {$item != ""} {
2640                        set val [format %-10s%9.3f \
2641                                [lindex $item 0],[lindex $item 1] \
2642                                [lindex $item 2]]
2643                        setexp $key $val $s 19
2644                    } else {
2645                        setexp $key " " $s 19
2646                    }
2647                }
2648                incr v0 3
2649            }
2650        }
2651        profile*-get {
2652            regsub profile $type {} term
2653            if {$term < 10} {
2654                set term " $term"
2655            }
2656            if {$number > 9} {
2657                set i [expr {($number+1)/10}]
2658            } else {
2659                set i " "
2660            }
2661            set key "LEQV PF$term  $i"
2662            # return nothing if no term exists
2663            if {![existsexp $key]} {return 0}
2664            # number of constraint lines
2665           
2666            set numline [string trim [string range [readexp ${key}] \
2667                    [expr {($number%10)*5}] [expr {4+(($number%10)*5)}]] ]
2668            if {$number == 0} {return $numline}
2669            set clist {}
2670            if {$number < 10} {
2671                set number " $number"
2672            }
2673            for {set i 1} {$i <= $numline} {incr i} {
2674                set key "LEQV PF${term}${number}$i"
2675                set line [readexp ${key}]
2676                foreach s {1 21 41} e {20 40 60} {
2677                    set seg [string range $line $s $e]
2678                    if {[string trim $seg] == ""} continue
2679                    # parse the string segment
2680                    set parse [regexp { *([0-9AL]+),([0-9AL]+) +([0-9.]+)} \
2681                            $seg junk phase hist mult]
2682                    # was parse successful
2683                    if {!$parse} {continue}
2684                    lappend clist [list $phase $hist $mult]
2685                }
2686            }
2687            return $clist
2688        }
2689        default {
2690            set msg "Unsupported constrinfo access: type=$type action=$action"
2691            tk_dialog .badexp "Error in readexp access" $msg error 0 OK
2692        }
2693
2694    }
2695}
2696
2697# read the default profile information for a histogram
2698# use: profdefinfo hist set# parm action
2699
2700#     proftype -- profile function number
2701#     profterms -- number of profile terms
2702#     pdamp -- damping value for the profile (*)
2703#     pcut -- cutoff value for the profile (*)
2704#     pterm$n -- profile term #n
2705#     pref$n -- refinement flag value for profile term #n (*)
2706
2707proc profdefinfo {hist set parm "action get"} {
2708    global expgui
2709    if {$hist < 10} {
2710        set key "HST  $hist"
2711    } else {
2712        set key "HST $hist"
2713    }
2714    switch -glob ${parm}-$action {
2715        proftype-get {
2716            set val [string range [readexp "${key}PRCF$set"] 0 4]
2717            if {$val == " "} {return 0}
2718            return $val
2719        }
2720        profterms-get {
2721            set val [string range [readexp "${key}PRCF$set"] 5 9]
2722            if {$val == " "} {return 0}
2723            return $val
2724        }
2725        pcut-get {
2726            return [string trim [string range [readexp "${key}PRCF$set"] 10 19]]
2727        }
2728        pdamp-get {
2729                set val [string range [readexp "${key}PRCF$set"] 24 24]
2730            if {$val == " "} {return 0}
2731            return $val
2732        }
2733        pterm*-get {
2734            regsub pterm $parm {} num
2735            set f1 [expr {15*(($num - 1) % 4)}]
2736            set f2 [expr {15*(1 + ($num - 1) % 4)-1}]
2737            set line  [expr {1 + ($num - 1) / 4}]
2738            return [string trim [string range [\
2739                        readexp "${key}PRCF${set}$line"] $f1 $f2] ]
2740        }
2741        pref*-get {
2742            regsub pref $parm {} num
2743            set f [expr {24+$num}]
2744            if {[string toupper [string range [readexp "${key}PRCF$set"] $f $f]] == "Y"} {
2745                return 1
2746            }
2747            return 0
2748        }
2749        default {
2750            set msg "Unsupported profdefinfo access: parm=$parm action=$action"
2751            tk_dialog .badexp "Code Error" $msg error 0 Exit
2752        }
2753    }
2754}
2755
2756# get March-Dollase preferred orientation information
2757# use MDprefinfo hist phase axis-number parm action value
2758#    ratio    -- ratio of xtallites in PO direction vs random (>1 for more)
2759#    fraction -- fraction in this direction, when more than one axis is used
2760#    h k & l  -- indices of P.O. axis
2761#    ratioref -- flag to vary ratio
2762#    fracref  -- flag to vary fraction
2763#    damp     -- damping value
2764#    type     -- model type (0 = P.O. _|_ to beam, 1 = || to beam)
2765#    new      -- creates a new record with default values (set only)
2766proc MDprefinfo {histlist phaselist axislist parm "action get" "value {}"} {
2767    foreach phase $phaselist hist $histlist axis $axislist {
2768        if {$phase == ""} {set phase [lindex $phaselist end]}
2769        if {$hist == ""} {set hist [lindex $histlist end]}
2770        if {$axis == ""} {set axis [lindex $axislist end]}
2771        if {$hist < 10} {
2772            set hist " $hist"
2773        }
2774        if {$axis > 9} {
2775            set axis "0"
2776        }
2777        set key "HAP${phase}${hist}PREFO${axis}"
2778        switch -glob ${parm}-$action {
2779            ratio-get {
2780                return [string trim [string range [readexp $key] 0 9]]
2781            }
2782            ratio-set {
2783                if ![validreal value 10 6] {return 0}
2784                setexp $key $value 1 10
2785            }
2786            fraction-get {
2787                return [string trim [string range [readexp $key] 10 19]]
2788            }
2789            fraction-set {
2790                if ![validreal value 10 6] {return 0}
2791                setexp $key $value 11 10
2792            }
2793            h-get {
2794                set h [string trim [string range [readexp $key] 20 29]]
2795                # why not allow negative h values?
2796                #               if {$h < 1} {return 0}
2797                return $h
2798            }
2799            h-set {
2800                if ![validreal value 10 2] {return 0}
2801                setexp $key $value 21 10
2802            }
2803            k-get {
2804                set k [string trim [string range [readexp $key] 30 39]]
2805                #               if {$k < 1} {return 0}
2806                return $k
2807            }
2808            k-set {
2809                if ![validreal value 10 2] {return 0}
2810                setexp $key $value 31 10
2811            }
2812            l-get {
2813                set l [string trim [string range [readexp $key] 40 49]]
2814                #if {$l < 1} {return 0}
2815                return $l
2816            }
2817            l-set {
2818                if ![validreal value 10 2] {return 0}
2819                setexp $key $value 41 10
2820            }
2821            ratioref-get {
2822                if {[string toupper \
2823                        [string range [readexp $key] 53 53]] == "Y"} {
2824                    return 1
2825                }
2826                return 0
2827            }
2828            ratioref-set {
2829                if $value {
2830                    setexp $key "Y" 54 1
2831                } else {
2832                    setexp $key "N" 54 1
2833                }
2834            }
2835            fracref-get {
2836                if {[string toupper \
2837                        [string range [readexp $key] 54 54]] == "Y"} {
2838                    return 1
2839                }
2840                return 0
2841            }
2842            fracref-set {
2843                if $value {
2844                    setexp $key "Y" 55 1
2845                } else {
2846                    setexp $key "N" 55 1
2847              }
2848            }
2849            damp-get {
2850                set val [string trim [string range [readexp $key] 59 59]]
2851                if {$val == " "} {return 0}
2852                return $val
2853            }
2854            damp-set {
2855                setexp $key $value 60 1
2856            }
2857            type-get {
2858                set val [string trim [string range [readexp $key] 64 64]]
2859                if {$val == " "} {return 0}
2860                return $val
2861            }
2862            type-set {
2863                # only valid settings are 0 & 1
2864                if {$value != "0" && $value != "1"} {set value "0"}
2865                setexp $key $value 65 1
2866            }
2867            new-set {
2868                makeexprec $key
2869                setexp $key \
2870                        {  1.000000  1.000000  0.000000  0.000000  1.000000   NN    0    0} \
2871                        1 68
2872            }
2873            default {
2874                set msg "Unsupported MDprefinfo access: parm=$parm action=$action"
2875                tk_dialog .badexp "Error in readexp" $msg error 0 Exit
2876            }
2877
2878        }
2879
2880    }
2881}
2882
2883# get list of defined atom types
2884proc AtmTypList {} {
2885    set natypes [readexp " EXPR  NATYP"]
2886    if {$natypes == ""} return
2887    set j 0
2888    set typelist {}
2889    for {set i 1} {$i <= $natypes} {incr i} {
2890        set key {this should never be matched}
2891        while {![existsexp $key]} {
2892            incr j
2893            if {$j > 99} {
2894                return $typelist
2895            } elseif {$j <10} {
2896                set key " EXPR ATYP $j"
2897            } else {
2898                set key " EXPR ATYP$j"
2899            }
2900        }
2901        lappend typelist [string trim [string range $::exparray($key) 2 9]]
2902    }
2903    return $typelist
2904}
2905
2906# read information about atom types
2907#     distrad    atomic distance search radius (get/set)
2908#     angrad     atomic angle search radius (get/set)
2909proc AtmTypInfo {parm atmtype "action get" "value {}"} {
2910    # first, search through the records to find the record matching the type
2911    set natypes [readexp " EXPR  NATYP"]
2912    if {$natypes == ""} return
2913    set j 0
2914    set typelist {}
2915    for {set i 1} {$i <= $natypes} {incr i} {
2916        set key {this should never be matched}
2917        while {![existsexp $key]} {
2918            incr j
2919            if {$j > 99} {
2920                return $typelist
2921            } elseif {$j <10} {
2922                set key " EXPR ATYP $j"
2923            } else {
2924                set key " EXPR ATYP$j"
2925            }
2926        }
2927        if {[string toupper $atmtype] == \
2928                [string toupper [string trim [string range $::exparray($key) 2 9]]]} break
2929        set key {}
2930    }
2931    if {$key == ""} {
2932        # atom type not found
2933        return {}
2934    }
2935    switch -glob ${parm}-$action {
2936        distrad-get {
2937            return [string trim [string range [readexp $key] 15 24]]
2938        }
2939        distrad-set {
2940            if ![validreal value 10 2] {return 0}
2941            setexp $key $value 16 10
2942        }
2943        angrad-get {
2944            return [string trim [string range [readexp $key] 25 34]]
2945        }
2946        angrad-set {
2947            if ![validreal value 10 2] {return 0}
2948            setexp $key $value 26 10
2949        }
2950        default {
2951            set msg "Unsupported AtmTypInfo access: parm=$parm action=$action"
2952            tk_dialog .badexp "Error in readexp" $msg error 0 Exit
2953        }
2954    }
2955}
2956# read default information about atom types (records copied to the .EXP file
2957# from the gsas/data/atomdata.dat file as AFAC ...
2958#     distrad returns a list of atom types (one or two letters) and
2959#                the corresponding distance
2960# note that these values are read only (no set option)
2961proc DefAtmTypInfo {parm} {
2962    set keys [array names ::exparray " AFAC *_SIZ"]
2963    set elmlist {}
2964    if {[llength $keys] <= 0} {return ""}
2965    foreach key $keys {
2966        lappend elmlist [string trim [string range $key 6 7]]
2967    }
2968    switch -glob ${parm} {
2969        distrad {
2970            set out {}
2971            foreach key $keys elm $elmlist {
2972                set val [string range $::exparray($key) 0 9]
2973                lappend out "$elm [string trim $val]"
2974            }
2975            return $out
2976        }
2977        angrad {
2978            set out {}
2979            foreach key $keys elm $elmlist {
2980                set val [string range $::exparray($key) 10 19]
2981                lappend out "$elm [string trim $val]"
2982            }
2983            return $out
2984        }
2985        default {
2986            set msg "Unsupported DefAtmTypInfo access: parm=$parm"
2987            tk_dialog .badexp "Error in readexp" $msg error 0 Exit
2988        }
2989    }
2990}
2991# write the .EXP file
2992proc expwrite {expfile} {
2993    global exparray
2994    set blankline \
2995     "                                                                        "
2996    set fp [open ${expfile} w]
2997    fconfigure $fp -translation crlf -encoding ascii
2998    set keylist [lsort [array names exparray]]
2999    # reorder the keys so that VERSION comes 1st
3000    set pos [lsearch -exact $keylist {     VERSION}]
3001    set keylist "{     VERSION} [lreplace $keylist $pos $pos]"
3002    foreach key $keylist {
3003        puts $fp [string range \
3004                "$key$exparray($key)$blankline" 0 79]
3005    }
3006    close $fp
3007}
3008
3009# history commands -- delete all but last $keep history records,
3010# renumber if $renumber is true
3011proc DeleteHistory {keep renumber} {
3012    global exparray
3013    foreach y [lrange [lsort -decreasing \
3014            [array names exparray {    HSTRY*}]] $keep end] {
3015        unset exparray($y)
3016    }
3017    if !$renumber return
3018    # renumber
3019    set i 0
3020    foreach y [lsort -increasing \
3021            [array names exparray {    HSTRY*}]] {
3022        set key [format "    HSTRY%3d" [incr i]]
3023        set exparray($key) $exparray($y)
3024        unset exparray($y)
3025    }
3026    # list all history
3027    #    foreach y [lsort -decreasing [array names exparray {    HSTRY*}]] {puts "$y $exparray($y)"}
3028}
3029
3030proc CountHistory {} {
3031    global exparray
3032    return [llength [array names exparray {    HSTRY*}]]
3033}
3034
3035# set the phase flags for histogram $hist to $plist
3036proc SetPhaseFlag {hist plist} {
3037    # make a 2 digit key -- hh
3038    if {$hist < 10} {
3039        set hh " $hist"
3040    } else {
3041        set hh $hist
3042    }
3043    set key "HST $hh NPHAS"
3044    set str {}
3045    foreach iph {1 2 3 4 5 6 7 8 9} {
3046        if {[lsearch $plist $iph] != -1} {
3047            append str {    1}
3048        } else {
3049            append str {    0}     
3050        }
3051    }
3052    setexp $key $str 1 68
3053}
3054
3055# erase atom $atom from phase $phase
3056# update the list of atom types, erasing the record if not needed.
3057proc EraseAtom {atom phase} {
3058    set type [atominfo $phase $atom type]
3059    if {$type == ""} return
3060    if {$atom < 10} {
3061        set key "CRS$phase  AT  $atom"
3062    } elseif {$atom < 100} {
3063        set key "CRS$phase  AT $atom"
3064    } else {
3065        set key "CRS$phase  AT$atom"
3066    }
3067    # delete the records for the atom
3068    global exparray
3069    foreach k [array names exparray ${key}*] {
3070        delexp $k
3071    }
3072    # change the number of atoms in the phase
3073    phaseinfo $phase natoms set [expr {[phaseinfo $phase natoms] -1}]
3074
3075    # now adjust numbers in "EXPR ATYP" records and delete, if needed.
3076    set natypes [readexp " EXPR  NATYP"]
3077    if {$natypes == ""} return
3078    set j 0
3079    for {set i 1} {$i <= $natypes} {incr i} {
3080        incr j
3081        if {$j <10} {
3082            set key " EXPR ATYP $j"
3083        } else {
3084            set key " EXPR ATYP$j"
3085        }
3086        while {![existsexp $key]} {
3087            incr j
3088            if {$j > 99} {
3089                return
3090            } elseif {$j <10} {
3091                set key " EXPR ATYP $j"
3092            } else {
3093                set key " EXPR ATYP$j"
3094            }
3095        }
3096        set keytype [string trim [string range $exparray($key) 2 9]]
3097        if {$type == $keytype} {
3098            # found the type record
3099            set val [string trim [string range $exparray($key) 10 14]]
3100            incr val -1
3101            # if this is the last reference, remove the record,
3102            # otherwise, decrement the counter
3103            if {$val <= 0} {
3104                incr natypes -1 
3105                validint natypes 5
3106                setexp " EXPR  NATYP" $natypes 1 5
3107                delexp $key
3108            } else {
3109                validint val 5
3110                setexp $key $val 11 5
3111            }
3112            return
3113        }
3114    }
3115}
3116
3117# compute equivalent anisotropic temperature factor for Uequiv
3118proc CalcAniso {phase Uequiv} {
3119    foreach var {a b c alpha beta gamma} {
3120        set $var [phaseinfo $phase $var]
3121    }
3122
3123    set G(1,1) [expr {$a * $a}]
3124    set G(2,2) [expr {$b * $b}]
3125    set G(3,3) [expr {$c * $c}]
3126    set G(1,2) [expr {$a * $b * cos($gamma*0.017453292519943)}]
3127    set G(2,1) $G(1,2)
3128    set G(1,3) [expr {$a * $c * cos($beta *0.017453292519943)}]
3129    set G(3,1) $G(1,3)
3130    set G(2,3) [expr {$b * $c * cos($alpha*0.017453292519943)}]
3131    set G(3,2) $G(2,3)
3132
3133    # Calculate the volume**2
3134    set v2 0.0
3135    foreach i {1 2 3} {
3136        set J [expr {($i%3) + 1}]
3137        set K [expr {(($i+1)%3) + 1}]
3138        set v2 [expr {$v2+ $G(1,$i)*($G(2,$J)*$G(3,$K)-$G(3,$J)*$G(2,$K))}]
3139    }
3140    if {$v2 > 0} {
3141        set v [expr {sqrt($v2)}]
3142        foreach i {1 2 3} {
3143            set i1 [expr {($i%3) + 1}]
3144            set i2 [expr {(($i+1)%3) + 1}]
3145            foreach j {1 2 3} {
3146                set j1 [expr {($j%3) + 1}]
3147                set j2 [expr {(($j+1)%3) + 1}]
3148                set C($j,$i) [expr {(\
3149                        $G($i1,$j1) * $G($i2,$j2) - \
3150                        $G($i1,$j2)  * $G($i2,$j1)\
3151                        )/ $v}]
3152            }
3153        }
3154        set A(1,2) [expr {0.5 * ($C(1,2)+$C(2,1)) / sqrt( $C(1,1)* $C(2,2) )}]
3155        set A(1,3) [expr {0.5 * ($C(1,3)+$C(3,1)) / sqrt( $C(1,1)* $C(3,3) )}]
3156        set A(2,3) [expr {0.5 * ($C(2,3)+$C(3,2)) / sqrt( $C(2,2)* $C(3,3) )}]
3157        foreach i {1 1 2} j {2 3 3} {
3158            set A($i,$j) [expr {0.5 * ($C($i,$j) + $C($j,$i)) / \
3159                    sqrt( $C($i,$i)* $C($j,$j) )}]
3160            # clean up roundoff
3161            if {abs($A($i,$j)) < 1e-5} {set A($i,$j) 0.0}
3162        }
3163    } else {
3164        set A(1,2) 0.0
3165        set A(1,3) 0.0
3166        set A(2,3) 0.0
3167    }
3168    return "$Uequiv $Uequiv $Uequiv \
3169            [expr {$Uequiv * $A(1,2)}] \
3170            [expr {$Uequiv * $A(1,3)}] \
3171            [expr {$Uequiv * $A(2,3)}]"
3172}
3173
3174# read/edit soft (distance) restraint info
3175#  parm:
3176#    weight -- histogram weight (factr) *
3177#    restraintlist -- list of restraints *
3178#  action: get (default) or set
3179#  value: used only with set
3180#  * =>  read+write supported
3181proc SoftConst {parm "action get" "value {}"} {
3182    set HST {}
3183    # look for RSN record
3184    set n 0
3185    for {set i 0} {$i < $::expmap(nhst)} {incr i} {
3186        set ihist [expr {$i + 1}]
3187        if {[expr {$i % 12}] == 0} {
3188            incr n
3189            set line [readexp " EXPR  HTYP$n"]
3190            if {$line == ""} {
3191                set msg "No HTYP$n entry for Histogram $ihist. This is an invalid .EXP file"
3192                tk_dialog .badexp "Error in readexp" $msg error 0 Exit
3193            }
3194            set j 0
3195        } else {
3196            incr j
3197        }
3198        if {[string range $line [expr 2+5*$j] [expr 5*($j+1)]] == "RSN "} {
3199            set HST $ihist
3200        }
3201    }
3202    if {$HST <=9} {
3203        set key "HST  $HST"
3204    } else {
3205        set key "HST $HST"
3206    }
3207    if {$HST == "" && $action == "set"} {
3208        # no RSN found need to add the soft constr. histogram
3209        # increment number of histograms
3210        set hst [string trim [string range [readexp { EXPR  NHST }] 0 4]]
3211        incr hst
3212        set HST $hst
3213        if ![validint hst 5] {return 0}
3214        setexp  { EXPR  NHST } $hst 1 5
3215        # add to EXPR HTYPx rec, creating if needed
3216        set n [expr { 1+ (($HST - 1) / 12) }]
3217        set key " EXPR  HTYP$n"
3218        if {[array names ::exparray $key] == ""} {
3219            makeexprec $key
3220        }
3221        setexp $key "RSN " [expr 3 + 5*(($HST-1) % 12)] 5
3222        # create other HST  xx recs
3223        if {$HST <=9} {
3224            set key "HST  $HST"
3225        } else {
3226            set key "HST $HST"
3227        }
3228        makeexprec "$key  HNAM"
3229        setexp "$key  HNAM" "Bond distance restraints" 3 24
3230        makeexprec "$key FACTR"
3231        makeexprec "$key NBNDS"
3232        mapexp
3233    } elseif {$HST == ""} {
3234        if $::expgui(debug) {puts "no restraints"}
3235        return "1"
3236    }
3237
3238    switch -glob ${parm}-$action {
3239        weight-get {
3240            return [string trim [string range [readexp "$key FACTR"] 0 14]]
3241        }
3242        weight-set {
3243            # update FACTR
3244            if ![validreal value 15 6] {return 0}
3245            setexp "$key FACTR" $value 1 15
3246        }
3247        restraintlist-get {
3248            set ncons [string trim [string range [readexp "$key NBNDS"] 0 4]]
3249            set conslist {}
3250            for {set i 1} {$i <= $ncons} {incr i} {
3251                set fi [string toupper [format %.4x $i]]
3252                set line [readexp "${key}BD$fi"]
3253                set const {}
3254                foreach len {3 5 5 3 3 3 3 3 6 6} {
3255                  set lenm1 [expr {$len - 1}]
3256                  lappend const [string trim [string range $line 0 $lenm1]]
3257                  set line [string range $line $len end]
3258                }
3259                lappend conslist $const
3260            }
3261            return $conslist
3262        }
3263        restraintlist-set {
3264            set num [llength $value]
3265            if ![validint num 5] {return 0}
3266            setexp "$key NBNDS" $num 1 5
3267            # delete all old records
3268            foreach i [array names ::exparray "${key}BD*"] {unset ::exparray($i)}
3269            set i 0
3270            foreach cons $value {
3271                incr i
3272                set fi [string toupper [format %.4x $i]]
3273                makeexprec "${key}BD$fi"
3274                set pos 1
3275                foreach num $cons len {3 5 5 3 3 3 3 3 -6 -6} {
3276                    if {$len > 0} {
3277                        validint num $len
3278                        setexp "${key}BD$fi" $num $pos $len
3279                    } else {
3280                        set len [expr abs($len)]
3281                        validreal num $len 3
3282                        setexp "${key}BD$fi" $num $pos $len
3283                    }
3284                    incr pos $len
3285                }
3286            }
3287        }
3288        default {
3289            set msg "Unsupported phaseinfo access: parm=$parm action=$action"
3290            tk_dialog .badexp "Error in readexp" $msg error 0 Exit
3291        }
3292    return 1
3293    }
3294}
3295
3296#======================================================================
3297# conversion routines
3298#======================================================================
3299
3300# convert x values to d-space
3301proc tod {xlist hst} {
3302    global expmap
3303    if {[string range $expmap(htype_$hst) 2 2] == "T"} {
3304        return [toftod $xlist $hst]
3305    } elseif {[string range $expmap(htype_$hst) 2 2] == "C"} {
3306        return [tttod $xlist $hst]
3307    } elseif {[string range $expmap(htype_$hst) 2 2] == "E"} {
3308        return [engtod $xlist $hst]
3309    } else {
3310        return {}
3311    }
3312}
3313
3314# convert tof to d-space
3315proc toftod {toflist hst} {
3316    set difc [expr {[histinfo $hst difc]/1000.}]
3317    set difc2 [expr {$difc*$difc}]
3318    set difa [expr {[histinfo $hst difa]/1000.}]
3319    set zero [expr {[histinfo $hst zero]/1000.}]
3320    set ans {}
3321    foreach tof $toflist {
3322        if {$tof == 0.} {
3323            lappend ans 0.
3324        } elseif {$tof == 1000.} {
3325            lappend ans 1000.
3326        } else {
3327            set td [expr {$tof-$zero}]
3328            lappend ans [expr {$td*($difc2+$difa*$td)/ \
3329                    ($difc2*$difc+2.0*$difa*$td)}]
3330        }
3331    }
3332    return $ans
3333}
3334
3335# convert two-theta to d-space
3336proc tttod {twotheta hst} {
3337    set lamo2 [expr {0.5 * [histinfo $hst lam1]}]
3338    set zero [expr [histinfo $hst zero]/100.]
3339    set ans {}
3340    set cnv [expr {acos(0.)/180.}]
3341    foreach tt $twotheta {
3342        if {$tt == 0.} {
3343            lappend ans 99999.
3344        } elseif {$tt == 1000.} {
3345            lappend ans 0.
3346        } else {
3347            lappend ans [expr {$lamo2 / sin($cnv*($tt-$zero))}]
3348        }
3349    }
3350    return $ans
3351}
3352
3353# convert energy (edx-ray) to d-space
3354# (note that this ignores the zero correction)
3355proc engtod {eng hst} {
3356    set lam [histinfo $hst lam1]
3357    set zero [histinfo $hst zero]
3358    set ans {}
3359    set v [expr {12.398/(2.0*[sind[expr ($lam/2.0)]])}]
3360    foreach e $eng {
3361        if {$e == 0.} {
3362            lappend ans 1000.
3363        } elseif {$e == 1000.} {
3364            lappend ans 0.
3365        } else {
3366            lappend ans [expr {$v/$e}]
3367        }
3368    }
3369    return $ans
3370}
3371
3372# convert x values to Q
3373proc toQ {xlist hst} {
3374    global expmap
3375    if {[string range $expmap(htype_$hst) 2 2] == "T"} {
3376        return [toftoQ $xlist $hst]
3377    } elseif {[string range $expmap(htype_$hst) 2 2] == "C"} {
3378        return [tttoQ $xlist $hst]
3379    } elseif {[string range $expmap(htype_$hst) 2 2] == "E"} {
3380        return [engtoQ $xlist $hst]
3381    } else {
3382        return {}
3383    }
3384}
3385# convert tof to Q
3386proc toftoQ {toflist hst} {
3387    set difc [expr {[histinfo $hst difc]/1000.}]
3388    set difc2 [expr {$difc*$difc}]
3389    set difa [expr {[histinfo $hst difa]/1000.}]
3390    set zero [expr {[histinfo $hst zero]/1000.}]
3391    set 2pi [expr {4.*acos(0.)}]
3392    set ans {}
3393    foreach tof $toflist {
3394        if {$tof == 0.} {
3395            lappend ans 99999.
3396        } elseif {$tof == 1000.} {
3397            lappend ans 0.
3398        } else {
3399            set td [expr {$tof-$zero}]
3400            lappend ans [expr {$2pi * \
3401                    ($difc2*$difc+2.0*$difa*$td)/($td*($difc2+$difa*$td))}]
3402        }
3403    }
3404    return $ans
3405}
3406
3407# convert two-theta to Q
3408proc tttoQ {twotheta hst} {
3409    set lamo2 [expr {0.5 * [histinfo $hst lam1]}]
3410    set zero [expr [histinfo $hst zero]/100.]
3411    set ans {}
3412    set cnv [expr {acos(0.)/180.}]
3413    set 2pi [expr {4.*acos(0.)}]
3414    foreach tt $twotheta {
3415        if {$tt == 0.} {
3416            lappend ans 0.
3417        } elseif {$tt == 1000.} {
3418            lappend ans 1000.
3419        } else {
3420            lappend ans [expr {$2pi * sin($cnv*($tt-$zero)) / $lamo2}]
3421        }
3422    }
3423    return $ans
3424}
3425# convert energy (edx-ray) to Q
3426# (note that this ignores the zero correction)
3427proc engtoQ {eng hst} {
3428    set lam [histinfo $hst lam1]
3429    set zero [histinfo $hst zero]
3430    set ans {}
3431    set v [expr {12.398/(2.0*[sind[expr ($lam/2.0)]])}]
3432    set 2pi [expr {4.*acos(0.)}]
3433    foreach e $eng {
3434        if {$e == 0.} {
3435            lappend ans 0.
3436        } elseif {$e == 1000.} {
3437            lappend ans 1000.
3438        } else {
3439            lappend ans [expr {$2pi * $e / $v}]
3440        }
3441    }
3442    return $ans
3443}
3444proc sind {angle} {
3445    return [expr {sin($angle*acos(0.)/90.)}]
3446}
3447
3448# convert d-space values to 2theta, TOF or KeV
3449proc fromd {dlist hst} {
3450    global expmap
3451    if {[string range $expmap(htype_$hst) 2 2] == "T"} {
3452        set difc [expr {[histinfo $hst difc]/1000.}]
3453        set difa [expr {[histinfo $hst difa]/1000.}]
3454        set zero [expr {[histinfo $hst zero]/1000.}]
3455        set ans {}
3456        foreach d $dlist {
3457            if {$d == 0.} {
3458                lappend ans 0.
3459            } elseif {$d == 1000.} {
3460                lappend ans 1000.
3461            } else {
3462                lappend ans [expr {$difc*$d + $difa*$d*$d + $zero}]
3463            }
3464        }
3465        return $ans
3466    } elseif {[string range $expmap(htype_$hst) 2 2] == "C"} {
3467        set lamo2 [expr {0.5 * [histinfo $hst lam1]}]
3468        set zero [expr [histinfo $hst zero]/100.]
3469        set ans {}
3470        set cnv [expr {180./acos(0.)}]
3471        foreach d $dlist {
3472            if {$d == 99999.} {
3473                lappend ans 0
3474            } elseif {$d == 0.} {
3475                lappend ans 1000.
3476            } else {
3477                lappend ans [expr {$cnv*asin($lamo2/$d) + $zero}]
3478            }
3479        }
3480        return $ans
3481    } elseif {[string range $expmap(htype_$hst) 2 2] == "E"} {
3482        set lam [histinfo $hst lam1]
3483        set zero [histinfo $hst zero]
3484        set v [expr {12.398/(2.0*[sind[expr ($lam/2.0)]])}]
3485        set ans {}
3486        foreach d $dlist {
3487            if {$d == 1000.} {
3488                lappend ans 0
3489            } elseif {$d == 0.} {
3490                lappend ans 1000.
3491            } else {
3492                lappend ans [expr {$v/$d}]
3493            }
3494        }
3495        return $ans
3496    } else {
3497        return {}
3498    }
3499}
3500
3501# convert Q values to 2theta, TOF or KeV
3502proc fromQ {Qlist hst} {
3503    global expmap
3504    if {[string range $expmap(htype_$hst) 2 2] == "T"} {
3505        set difc [expr {[histinfo $hst difc]/1000.}]
3506        set difa [expr {[histinfo $hst difa]/1000.}]
3507        set zero [expr {[histinfo $hst zero]/1000.}]
3508        set ans {}
3509        foreach Q $Qlist {
3510            if {$Q == 0.} {
3511                lappend ans 1000.
3512            } elseif {$Q == 99999.} {
3513                lappend ans 1000.
3514            } else {
3515                set d [expr {4.*acos(0.)/$Q}]
3516                lappend ans [expr {$difc*$d + $difa*$d*$d + $zero}]
3517            }
3518        }
3519        return $ans
3520    } elseif {[string range $expmap(htype_$hst) 2 2] == "C"} {
3521        set lamo4pi [expr {[histinfo $hst lam1]/(8.*acos(0.))}]
3522        set zero [expr [histinfo $hst zero]/100.]
3523        set ans {}
3524        set cnv [expr {180./acos(0.)}]
3525        foreach Q $Qlist {
3526            if {$Q == 0.} {
3527                lappend ans 0
3528            } elseif {$Q == 1000.} {
3529                lappend ans 1000.
3530            } else {
3531                lappend ans [expr {$cnv*asin($Q*$lamo4pi) + $zero}]
3532            }
3533        }
3534        return $ans
3535    } elseif {[string range $expmap(htype_$hst) 2 2] == "E"} {
3536        set lam [histinfo $hst lam1]
3537        set zero [histinfo $hst zero]
3538        set v [expr {12.398/(2.0*[sind[expr ($lam/2.0)]])}]
3539        set ans {}
3540        set 2pi [expr {4.*acos(0.)}]
3541        foreach Q $Qlist {
3542            if {$Q == 1000.} {
3543                lappend ans 0
3544            } elseif {$Q == 0.} {
3545                lappend ans 1000.
3546            } else {
3547                lappend ans [expr {$Q * $v/$2pi}]
3548            }
3549        }
3550        return $ans
3551    } else {
3552        return {}
3553    }
3554}
3555#============================================================================
3556# rigid body EXP editing routines (to move into readexp.tcl)
3557# RigidBodyList -- returns a list of the defined rigid body types
3558# SetRigidBodyVar -- set variables and damping for rigid body type multipliers
3559# ReadRigidBody  -- # of times a body is mapped, scaling factors, var #s & coordinates
3560# RigidBodyMappingList - return a list instances where a RB is mapped in phase
3561# RigidBodyEnableTLS -- Enable or Disable TLS use for a rigid body mapping
3562# RigidBodySetTLS  -- change the TLS values for a rigid body mapping
3563# RigidBodySetDamp -- change the damping values for a rigid body mapping
3564# RigidBodyVary    -- set refinement variable numbers for a rigid body mapping
3565# RigidBodyTLSVary -- set TLS refinement variable nums for a rigid body mapping
3566# AddRigidBody -- defines a new rigid body type
3567# DeleteRigidBody -- remove a rigid body definition
3568# ReplaceRigidBody -- replaces a previous rigid body type
3569# ReadRigidBodyMapping  -- get parameters for a rigid body mapping
3570# MapRigidBody -- map a rigid body type into a phase
3571# EditRigidBodyMapping -- change the parameters in a rigid body mapping
3572# UnMapRigidBody --remove a rigid body constraint by removing a RB "instance"
3573#----- note that these older routines should not be used ------
3574# RigidBodyCount -- returns the number of defined rigid bodies (body types)
3575#    use RigidBodyList instead
3576# RigidBodyMappingCount -- # of times a rigid body is mapped in phase
3577#    use RigidBodyMappingList instead
3578#============================================================================
3579# returns the number of defined rigid bodies
3580proc RigidBodyCount {} {
3581    set n [string trim [readexp "RGBD  NRBDS"]]
3582    if {$n == ""} {
3583        set n 0
3584    }
3585    return $n
3586}
3587
3588# returns a list of the defined rigid body types
3589proc RigidBodyList {} {
3590    set n [string trim [readexp "RGBD  NRBDS"]]
3591    if {$n == ""} {
3592        set n 0
3593    }
3594    set rblist {}
3595    foreach rbnum {1 2 3 4 5 6 7 8 9 10 11 12 13 14 15} {
3596        set value $rbnum
3597        validint value 2
3598        set key "RGBD${value}"
3599        if {[existsexp "$key NATR "]} {
3600            lappend rblist $rbnum
3601        }
3602        if {[llength $rblist] == $n} break
3603    }
3604    return $rblist
3605}
3606
3607# ReadRigidBody provides all information associated with a rigid body type
3608#  rbnum is the rigid body type number
3609# it returns two items:
3610#   the number of times the rigid body is mapped
3611#   a list containing an element for each scaling factor in rigid body #rbnum.
3612# in each element there are four items:
3613#    a multiplier value for the rigid body coordinates
3614#    a damping value (0-9) for the refinement of the multiplier
3615#    a variable number if the multiplier will be refined
3616#    a list of cartesian coordinates coordinates
3617# each cartesian coordinate contains 4 items: x,y,z and a label
3618#  note that the label is present only when the RB is created in EXPGUI and is
3619#  not used in GSAS.
3620proc ReadRigidBody {rbnum} {
3621    if {[lsearch [RigidBodyList] $rbnum] == -1} {
3622        return ""
3623    }
3624    set value $rbnum
3625    validint value 2
3626    set key "RGBD${value}"
3627    set n [string trim [string range [readexp "$key NATR"] 0 4]]
3628    set used [string trim [string range [readexp "$key NBDS"] 0 4]]
3629    set nmult [string trim [string range [readexp "$key NSMP"] 0 4]]
3630    set out {}
3631    for {set i 1} {$i <= $nmult} {incr i} {
3632        set line [readexp "${key}${i}PARM"]
3633        set mult [string trim [string range $line 0 9]]
3634        set var [string trim [string range $line 10 14]]
3635        set damp [string trim [string range $line 15 19]]
3636        set coordlist {}
3637        for {set j 1} {$j <= $n} {incr j} {
3638            set value $j
3639            validint value 3
3640            set line [readexp "${key}${i}SC$value"]
3641            set x [string trim [string range $line 0 9]]
3642            set y [string trim [string range $line 10 19]]
3643            set z [string trim [string range $line 20 29]]
3644            set lbl [string trim [string range $line 30 39]]
3645            lappend coordlist [list $x $y $z $lbl]
3646        }
3647        lappend out [list $mult $damp $var $coordlist]
3648    }
3649    return [list $used $out]
3650}
3651
3652# SetRigidBodyVar
3653#   rbnum is the rigid body type number
3654#   varnumlist is a list of variable numbers
3655#      note that if this list is shorter than the number of actual multipliers
3656#      for the body, the unspecified variable will not be changed
3657#   damplist   is a list of damping values (0-9)
3658#      note that if the damplist is shorter than the number of actual multipliers
3659#      the unspecified values are not changed
3660#  SetRigidBodVar 2 {1 2 3} {}
3661#       will vary the (first 3) translations in body #3 and will not change the
3662#       damping values
3663#  SetRigidBodVar 3 {} {0 0 0}
3664#       will not change variable settings but will change the (first 3) damping values
3665#  SetRigidBodVar 4 {11 11} {8 8}
3666#      changes both variable numbers and damping at the same time
3667# Nothing is returned
3668proc SetRigidBodyVar {rbnum varnumlist damplist} {
3669    if {[lsearch [RigidBodyList] $rbnum] == -1} {
3670        return ""
3671    }
3672    set value $rbnum
3673    validint value 2
3674    set key "RGBD${value}"
3675    set nmult [string trim [string range [readexp "$key NSMP"] 0 4]]
3676    for {set i 1} {$i <= $nmult} {incr i} {
3677        set j $i
3678        incr j -1
3679        set var [lindex $varnumlist $j]
3680        if {$var != ""} {
3681            validint var 5
3682            setexp "${key}${i}PARM" $var 11 15
3683        }
3684        set damp [lindex $damplist $j]
3685        if {$damp != ""} {
3686            if {$damp > 9} {set damp 9}
3687            if {$damp < 0} {set damp 0}
3688            validint damp 5
3689        }
3690        setexp "${key}${i}PARM" $damp 16 20
3691    }
3692}
3693
3694
3695# return the number of times rigid body $bodytyp is mapped in phase $phase
3696proc RigidBodyMappingCount {phase bodytyp} {
3697    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1]"
3698    if {! [existsexp "$key  NBDS"]} {return 0}
3699    set n [string trim [readexp "$key  NBDS"]]
3700    if {$n == ""} {
3701        set n 0
3702    }
3703    return $n
3704}
3705# return a list of the instances where rigid body $bodytyp is mapped in phase $phase
3706proc RigidBodyMappingList {phase bodytyp} {
3707    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1]"
3708    if {! [existsexp "$key  NBDS"]} {return {}}
3709    set n [string trim [readexp "$key  NBDS"]]
3710    if {$n == ""} {
3711        set n 0
3712    }
3713    set rblist {}
3714    foreach rbnum {1 2 3 4 5 6 7 8 9 10 11 12 13 14 15} {
3715        set value $rbnum
3716        validint value 2
3717        set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $rbnum 1]"
3718        if {[existsexp "$key  NDA"]} {
3719            lappend rblist $rbnum
3720        }
3721        if {[llength $rblist] == $n} break
3722    }
3723    return $rblist
3724}
3725
3726
3727
3728# reads rigid body mapping parameters for phase ($phase), body type # ($bodytyp) and instance # ($num)
3729# returns a list of items (most lists) as follows:
3730#   1) sequence # of first atom in body
3731#   2) origin of body in fractional coordinates (3 elements)
3732#   3) Euler angles as 6 pairs of numbers (see below)
3733#   4) variable numbers for the 9 position variables (origin followed by rotations)
3734#   5) damping vals for the 9 position variables (origin followed by rotations)
3735#   6) the TLS values, in order below (empty list if TLS is not in use)
3736#   7) the variable numbers for each TLS values, in order below (or empty)
3737#   8) three damping values for the T, L and S terms.
3738# returns an empty list if no such body exists.
3739#
3740# Euler angles are a list of axes and angles to rotate:
3741#   { {axis1 angle1} {axis2 angle2} ...}
3742# where axis1,... can be 1, 2 or 3 corresponding to the cartesian X, Y or Z axes
3743#
3744# The 20 TLS terms are ordered:
3745#    T11, T22, T33, T12, T13, T23
3746#    L11, L22, L33, L12, L13, L23
3747#    S12, S13, S21, S23, S31, S32, SAA, SBB
3748#
3749proc ReadRigidBodyMapping {phase bodytyp num} {
3750    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
3751        return ""
3752    }
3753    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
3754    set first [string trim [string range [readexp "$key  NDA"] 0 4]]
3755    set line [readexp "$key BDFL"]
3756    set varlist {}
3757    set damplist {}
3758    foreach i {0 1 2 3 4 5 6 7 8} {
3759        lappend varlist [string trim [string range $line [expr {5*$i}] [expr {4 + 5*$i}] ]]
3760        lappend damplist [string trim [string range $line [expr {45 + $i}] [expr {45 + $i}] ]]
3761    }
3762    set TLSdamplist {}
3763    foreach i {54 55 56} {
3764        lappend TLSdamplist [string trim [string range $line $i $i ]]
3765    }
3766    set line [readexp "${key} BDLC"]
3767    set x [string trim [string range $line 0 9]]
3768    set y [string trim [string range $line 10 19]]
3769    set z [string trim [string range $line 20 29]]
3770    set origin [list $x $y $z]
3771    set line [readexp "${key} BDOR"]
3772    set rotations {}
3773    foreach i {0 10 20 30 40 50} {
3774        set angle [string trim [string range $line $i [expr {$i+7}]]]
3775        set axis [string trim [string range $line [expr {$i+8}] [expr {$i+9}]]]
3776        lappend rotations [list $angle $axis]
3777    }
3778    set TLS [string trim [string range [readexp "${key} LSTF"] 0 4]]
3779    set tlsvars {}
3780    set tlsvals {}
3781    if {$TLS != 0} {
3782        set line [readexp "${key}TLSF1"]
3783        for {set j 0} {$j < 20} {incr j} {
3784            set var [string trim [string range $line [expr {3*$j}] [expr {3*$j+2}]]]
3785            if {$var == ""} {set var 0}
3786            lappend tlsvars $var
3787        }
3788        for {set j 0} {$j < 20} {incr j} {
3789            set i 0
3790            if {$j == 0} {
3791                set i 1
3792            } elseif {$j == 8} {
3793                set i 2
3794            } elseif {$j == 16} {
3795                set i 3
3796            }
3797            if {$i != 0} {
3798                set line [readexp "${key}TLSP$i"]
3799                set i 0
3800                set j1 0
3801                set j2 7
3802            } else {
3803                incr j1 8
3804                incr j2 8
3805            }
3806            set val [string trim [string range $line $j1 $j2]]
3807            if {$val == ""} {set val 0}
3808            lappend tlsvals $val
3809        }
3810    }
3811    return [list $first $origin $rotations $varlist $damplist $tlsvals $tlsvars $TLSdamplist]
3812}
3813
3814# Control TLS representation for phase, body # and instance number of a Rigid body mapping
3815#   for mapping with phase ($phase), body type # ($bodytyp) and instance # ($num)
3816# Enable TLS use if TLS is non-zero (true). Disable if zero
3817proc RigidBodyEnableTLS {phase bodytyp num TLS} {
3818    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
3819        return ""
3820    }
3821    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
3822    if {$TLS} {
3823        setexp "${key} LSTF" [format "%5d" 1] 1 5
3824        if {![existsexp "${key}TLSF1"]} {makeexprec "${key}TLSF1"}
3825        if {![existsexp "${key}TLSP1"]} {
3826            makeexprec "${key}TLSP1"
3827            set str {}
3828            foreach v {.01 .01 .01 0 0 0 0 0} d {4 4 4 4 4 4 2 2} {
3829                validreal v 8 $d
3830                append str $v
3831            }
3832            setexp "${key}TLSP1" $str 1 64
3833        }
3834        if {![existsexp "${key}TLSP2"]} {
3835            makeexprec "${key}TLSP2"
3836            set str {}
3837            set v 0
3838            foreach d {2 2 2 2 4 4 4 4} {
3839                validreal v 8 $d
3840                append str $v
3841            }
3842            setexp "${key}TLSP2" $str 1 64
3843        }
3844        if {![existsexp "${key}TLSP3"]} {
3845            makeexprec "${key}TLSP3"
3846            set str {}
3847            set v 0
3848            foreach d {4 4 4 4} {
3849                validreal v 8 $d
3850                append str $v
3851            }
3852            setexp "${key}TLSP3" $str 1 64
3853        }
3854    } else {
3855        setexp "${key} LSTF" [format "%5d" 0] 1 5
3856    }
3857    return 1
3858}
3859
3860# Control the TLS values for Rigid body mapping for mapping with
3861#    phase ($phase), body type # ($bodytyp) and instance # ($num)
3862# set the 20 TLS values to the values in TLSvals
3863# There must be exactly 20 TLS terms, which are ordered:
3864#    T11, T22, T33, T12, T13, T23
3865#    L11, L22, L33, L12, L13, L23
3866#    S12, S13, S21, S23, S31, S32, SAA, SBB
3867proc RigidBodySetTLS {phase bodytyp num TLSvals} {
3868    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
3869        return ""
3870    }
3871    if {[llength $TLSvals] != 20} {return ""}
3872    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
3873    set TLS [string trim [string range [readexp "${key} LSTF"] 0 4]]
3874    if {$TLS == 0} {return ""}
3875    if {![existsexp "${key}TLSF1"]} {makeexprec "${key}TLSF1"}
3876    foreach n {1 2 3} {
3877        if {![existsexp "${key}TLSP$n"]} {makeexprec "${key}TLSP$n"}
3878    }
3879    set str {}
3880    set n 1
3881    set i 0
3882    foreach v $TLSvals d {4 4 4 4 4 4 2 2 2 2 2 2 4 4 4 4 4 4 4 4} {
3883        incr i
3884        validreal v 8 $d
3885        append str $v
3886        if {$i == 8} {
3887            set i 0
3888            setexp "${key}TLSP$n" $str 1 64
3889            incr n
3890            set str {}
3891        }
3892    }
3893    setexp "${key}TLSP$n" $str 1 64
3894    return 1
3895}
3896
3897# set damping values for a Rigid body mapping
3898#   for mapping with phase ($phase), body type # ($bodytyp) and instance # ($num)
3899# there must be 9 damping values in RBdamp for the 9 position variables (origin followed by rotations)
3900# Use of TLSdamp is optional, but to be used, TLS representation must be enabled and there must be
3901# three damping terms (for all T terms; for all L terms and for all S terms)
3902proc RigidBodySetDamp {phase bodytyp num RBdamp "TLSdamp {}"} {
3903    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
3904        return ""
3905    }
3906    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
3907    if {[llength $RBdamp] != 9} {return ""}
3908    set str {}
3909    foreach v $RBdamp {
3910        if {[validint v 1] != 1} {set v " "}
3911        append str $v
3912    }
3913    setexp "$key BDFL" $str 46 9
3914    set TLS [string trim [string range [readexp "${key} LSTF"] 0 4]]
3915    if {$TLS != 0 &&  [llength $TLSdamp] == 3} {
3916        set str {}
3917        foreach v $TLSdamp {
3918        if {[validint v 1] != 1} {set v " "}
3919            append str $v
3920        }
3921        setexp "$key BDFL" $str 55 3
3922    }
3923    return 1
3924}
3925
3926# set refinement variable numbers for a Rigid body mapping
3927#   for mapping with phase ($phase), body type # ($bodytyp) and instance # ($num)
3928# there must be 9 variable values in RBvar for the 9 position variables (origin followed by rotations)
3929# note that the variable values should be unique integers
3930proc RigidBodyVary {phase bodytyp num RBvar} {
3931    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
3932        return ""
3933    }
3934    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
3935    if {[llength $RBvar] != 9} {return ""}
3936    set str {}
3937    foreach v $RBvar {
3938        if {[validint v 5] != 1} {set v " "}
3939        append str $v
3940    }
3941    setexp "$key BDFL" $str 1 45   
3942}
3943
3944# set TLS refinement variable numbers for a Rigid body mapping
3945#   for mapping with phase ($phase), body type # ($bodytyp) and instance # ($num)
3946# there must be 20 variable values in TLSvar for the 20 parameters:
3947#    T11, T22, T33, T12, T13, T23
3948#    L11, L22, L33, L12, L13, L23
3949#    S12, S13, S21, S23, S31, S32, SAA, SBB
3950# note that the variable values should be unique integers
3951proc RigidBodyTLSVary {phase bodytyp num TLSvar} {
3952    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $num] == -1} {
3953        return ""
3954    }
3955    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $num 1]"
3956    if {[llength $TLSvar] != 20} {return ""}
3957    set TLS [string trim [string range [readexp "${key} LSTF"] 0 4]]
3958    if {$TLS == 0} {return ""}
3959    set str {}
3960    foreach v $TLSvar {
3961        if {[validint v 3] != 1} {set v " "}
3962        append str $v
3963    }
3964    setexp "${key}TLSF1" $str 1 60
3965
3966# AddRigidBody: add a new rigid body definition into the .EXP file
3967# arguments are:
3968#   multlist: defines a list of multipliers for each set of coordinates. In the
3969#             simplest case this will be {1}
3970#   coordlist: a nested list of coordinates such as { { {0 0 0} {.1 .1 .1} {.2 .2 .2} } }
3971# note that when the length of multlist > 1 then coordlist must have the same length.
3972# for input where
3973#     multlist = {s1 s2} and
3974#     coordlist = { { {0 0 0} {1 1 0} {.0 .0 .0} ...}
3975#                     {0 0 0} {1 1 0} {2 1 2} ...}
3976#                 }
3977# the cartesian coordinates are defined from the input as
3978#    atom 1 = s1 * (0,0,0) + s2*(0,0,0) [ = (0,0,0)]
3979#    atom 2 = s1 * (1,1,0) + s2*(1,1,0) [ = (s1+s2) * (1,1,0)]
3980#    atom 3 = s1 * (0,0,0) + s2*(2,1,2) [ = s2 * (2,1,2)]
3981#    ...
3982# Returns the number of the rigid body that has been created
3983proc AddRigidBody {multlist coordlist} {
3984    # find the first unused body #
3985    foreach rbnum {1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16} {
3986        set value $rbnum
3987        validint value 2
3988        set key "RGBD${value}"
3989        if {! [existsexp "$key NATR "]} {break}
3990    }
3991    # did we go too far?
3992    if {$rbnum == 16} {return ""}
3993    # increment the RB counter
3994    set n [string trim [readexp "RGBD  NRBDS"]]
3995    if {$n == ""} {
3996        makeexprec "RGBD  NRBDS"
3997        set n 0
3998    }
3999    incr n
4000    validint n 5
4001    setexp "RGBD  NRBDS" $n 1 5
4002    SetRigidBody $rbnum $multlist $coordlist
4003    return $rbnum
4004}
4005
4006# DeleteRigidBody: remove a rigid body definition from the .EXP file
4007# The body may not be mapped. I am not sure if GSAS allows more than 9 bodies,
4008# but if it does, the simplifed approach used here will fail, so this
4009# is not allowed.
4010# Input:
4011#   Rigid body number
4012# Returns:
4013#   1 on success
4014#   -1 if the body number is 11 or greater
4015#   -2 if the body is mapped
4016#   -3 if the body is not defined
4017proc DeleteRigidBody {rbnum} {
4018    # can't delete bodies with numbers higher than 10, since the key prefix
4019    # (RGBD11... will overlap with rigid body instance records, which would be
4020    # deleted below
4021    if {$rbnum > 10} {
4022        return -1
4023    }
4024    set value $rbnum
4025    validint value 2
4026    set key "RGBD${value}"
4027    if {![existsexp "$key NATR "]} {
4028        return -2
4029    }
4030    # make sure the body is not mapped
4031    if {[string trim [string range [readexp "$key NBDS"] 0 4]] != 0} {
4032        return -3
4033    }
4034    # delete the records starting with "RGBD x" or "RGBD10"
4035    foreach key [array names ::exparray "${key}*"] {
4036        #puts $key
4037        delexp $key
4038    }
4039    # decrement the RB counter
4040    set n [string trim [readexp "RGBD  NRBDS"]]
4041    if {$n == ""} {
4042        set n 0
4043    }
4044    incr n -1
4045    validint n 5
4046    if {$n > 0} {
4047        setexp "RGBD  NRBDS" $n 1 5
4048    } else {
4049        delexp "RGBD  NRBDS"
4050    }
4051    return 1
4052}
4053
4054# ReplaceRigidBody: replace all the information for rigid body #rbnum
4055# Works the sames as AddRigidBody (see above) except that the rigid body is replaced rather
4056# than added.
4057# Note that count of the # of times the body is used is preserved
4058proc ReplaceRigidBody {rbnum multlist coordlist {varlist ""} {damplist ""}} {
4059    set value $rbnum
4060    validint value 2
4061    set key "RGBD${value}"
4062    set line [readexp "$key NBDS"]
4063    foreach key [array names ::exparray "${key}*"] {
4064        #puts $key
4065        delexp $key
4066    }
4067    SetRigidBody $rbnum $multlist $coordlist $varlist $damplist
4068    setexp "$key NBDS" $line 1 68
4069}
4070
4071# Edit the parameters for rigid body #rbnum
4072# (normally called from ReplaceRigidBody or AddRigidBody)
4073proc SetRigidBody {rbnum multlist coordlist {varlist ""} {damplist ""}} {
4074    set value $rbnum
4075    validint value 2
4076    set key "RGBD${value}"
4077    # number of atoms
4078    set value [llength [lindex $coordlist 0]]
4079    validint value 5
4080    makeexprec "$key NATR"
4081    setexp "$key NATR" $value 1 5
4082    # number of times used
4083    set value 0
4084    validint value 5
4085    makeexprec "$key NBDS"
4086    setexp "$key NBDS" $value 1 5
4087    # number of coordinate matrices
4088    set value [llength $multlist]
4089    validint value 5
4090    makeexprec "$key NSMP"
4091    setexp "$key NSMP" $value 1 5
4092    set i 0
4093    foreach mult $multlist coords $coordlist {
4094        set var [lindex $varlist $i]
4095        if {$var == ""} {set var 0}
4096        set damp [lindex $damplist $i]
4097        if {$damp == ""} {set damp 0}
4098        incr i
4099        makeexprec "${key}${i}PARM"
4100        setexp "${key}${i}PARM" [format "%10.5f%5d%5d" $mult $var $damp] 1 20
4101        set j 0
4102        foreach item $coords {
4103            #puts $item
4104            incr j
4105            set value $j
4106            validint value 3
4107            makeexprec "${key}${i}SC$value"
4108            if {[llength $item] == 4} {
4109                setexp "${key}${i}SC$value" [eval format "%10.6f%10.6f%10.6f%10s" $item] 1 40
4110            } elseif {[llength $item] == 3} {
4111                setexp "${key}${i}SC$value" [eval format "%10.6f%10.6f%10.6f" $item] 1 30
4112            } else {
4113                return -code 3 "Invalid number of coordinates"
4114            }
4115        }
4116    }
4117}
4118
4119# convert a decimal to the GSAS hex encoding with a field $digits long.
4120proc ToHex {num digits} {
4121    return [string toupper [format "%${digits}x" $num]]
4122}
4123
4124# convert a GSAS hex encoding to a decimal integer
4125proc FromHex {hex} {
4126    return [scan $hex "%x"]
4127}
4128
4129# MapRigidBody: define an "instance" of a rigid body: meaning that the coordinates
4130# (and optionally U values) for a set of atoms will be generated from the rigid body
4131# arguments:
4132#   phase: phase number (1-9)
4133#   bodytyp: number of rigid body (1-15) as returned from AddRigidBody
4134#   firstatom: sequence number of the first atom in phase (note that atoms may
4135#              not be numbered sequentially)
4136#   position: list of three fractional coordinates for the origin of the rigid body coordinates
4137#   angles: list of 3 angles to rotate the rigid body coordinates around x, y, z of the
4138#           cartesian system before the body is translated to position.
4139# returns the instance # (number of times body $bodytyp has been used in phase $phase)
4140proc MapRigidBody {phase bodytyp firstatom position angles} {
4141    # find the first unused body # for this phase & type
4142    foreach rbnum {1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16} {
4143        set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $rbnum 1]"
4144        if {! [existsexp "$key  NDA"]} {break}
4145    }
4146    # did we go too far?
4147    if {$rbnum == 16} {return ""}
4148    # increment number of mapped bodies of this type overall
4149    set value $bodytyp
4150    validint value 2
4151    set key "RGBD${value}"
4152    set used [string trim [string range [readexp "$key NBDS"] 0 4]]
4153    incr used
4154    set value $used
4155    validint value 5
4156    setexp "$key NBDS" $value 1 5
4157    # increment number of mapped bodies of this type in this phase
4158    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1]"
4159    if {[existsexp "$key  NBDS"]} {
4160        set used [string trim [string range [readexp "$key  NBDS"] 0 4]]
4161    } else {
4162        makeexprec "$key  NBDS"
4163        set used 0
4164    }
4165    incr used
4166    set value $used
4167    validint value 5
4168    setexp "$key  NBDS" $value 1 5
4169    # now write the mapping parameters
4170    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $rbnum 1]"
4171    set value $firstatom
4172    validint value 5
4173    makeexprec "$key  NDA"
4174    setexp "$key  NDA" $value 1 5
4175    set l1 {}
4176    set l2 {}
4177    for {set i 0} {$i < 9} {incr i} {
4178        append l1 [format %5d 0]
4179        append l2 [format %1d 0]
4180    }
4181    makeexprec "$key BDFL"
4182    setexp "$key BDFL" $l1$l2 1 54
4183    makeexprec "${key} BDLC"
4184    setexp "${key} BDLC" [eval format "%10.6f%10.6f%10.6f" $position] 1 30
4185    makeexprec "${key} BDOR"
4186    set l1 {}
4187    foreach val "$angles 0 0 0" dir "1 2 3 1 1 1" {
4188        append l1 [format "%8.2f%2d" $val $dir]
4189    }
4190    setexp "${key} BDOR" $l1 1 60
4191    makeexprec "${key} LSTF"
4192    setexp "${key} LSTF" [format "%5d" 0] 1 5
4193    return $rbnum
4194}
4195
4196# EditRigidBodyMapping: edit parameters that define an "instance" of a rigid body (see MapRigidBody)
4197# arguments:
4198#   phase: phase number (1-9)
4199#   bodytyp: number of rigid body (1-15) as returned from AddRigidBody
4200#   bodynum: instance number, as returned by MapRigidBody
4201#   position: list of three fractional coordinates for the origin of the rigid body coordinates
4202#   angles: list of 3 angles to rotate the rigid body coordinates around x, y, z of the
4203#           cartesian system before the body is translated to position.
4204#
4205proc EditRigidBodyMapping {phase bodytyp bodynum position angles} {
4206    # number of bodies of this type in this phase
4207    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $bodynum 1]"
4208    setexp "${key} BDLC" [eval format "%10.6f%10.6f%10.6f" $position] 1 30
4209    set l1 {}
4210    foreach val "$angles 0 0 0" dir "1 2 3 1 1 1" {
4211        append l1 [format "%8.2f%2d" $val $dir]
4212    }
4213    setexp "${key} BDOR" $l1 1 60
4214}
4215
4216# UnMapRigidBody: remove a rigid body constraint by removing a RB "instance"
4217# (undoes MapRigidBody)
4218# arguments:
4219#   phase: phase number (1-9)
4220#   bodytyp: number of rigid body (1-15) as returned from AddRigidBody
4221#   bodynum: instance number, as returned by MapRigidBody
4222proc UnMapRigidBody {phase bodytyp mapnum} {
4223    if {[lsearch [RigidBodyMappingList $phase $bodytyp] $mapnum] == -1} {
4224        return ""
4225    }
4226    # decrement number of mapped bodies of this type overall
4227    set value $bodytyp
4228    validint value 2
4229    set key "RGBD${value}"
4230    set used [string trim [string range [readexp "$key NBDS"] 0 4]]
4231    incr used -1
4232    set value $used
4233    validint value 5
4234    setexp "$key NBDS" $value 1 5
4235    # decrement number of mapped bodies of this type in this phase
4236    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1]"
4237    if {[existsexp "$key  NBDS"]} {
4238        set used [string trim [string range [readexp "$key  NBDS"] 0 4]]
4239    } else {
4240        set used 0
4241    }
4242    incr used -1
4243    set value $used
4244    validint value 5
4245    if {$used > 0} {
4246        setexp "$key  NBDS" $value 1 5
4247    } else {
4248        delexp "$key  NBDS"
4249    }
4250    # now delete the mapping parameter records
4251    set key "RGBD[ToHex $phase 1][ToHex $bodytyp 1][ToHex $mapnum 1]"
4252    foreach key [array names ::exparray "${key}*"] {
4253        delexp $key
4254    }
4255    return $used
4256}
4257
Note: See TracBrowser for help on using the repository browser.