Liberty BASIC Community Forum
Liberty BASIC Programming Discussions >> Novice >> simple math but confusing coding

simple math but confusing coding
Post by blessed2b3rd on Oct 20th, 2015, 12:08pm

Hi everyone,

I know how you guys like to use coding for real life issues, so I have one for work that I'm just not getting quite right. The math I can do manually its the coding loops and logic thats catching me out.

Context: In any mining company there will be a drillhole database that is used to store records of what the geologist has logged or what values the laboratory has returned. It is these values that are used to calculate how much ore is in the ground. A drillhole can be anything from 1 meter to a kilometer or even more. The deeper the ore deposit the deeper the holes will go. A geologist's job is to look at the rock samples in the core-box and describe what kind of rock it is and its interval width. This is known as logging. The geologist will 'log' all of the intervals with a description in either some software or excel. Usually giving each interval a rock code (the most common codes from a drop-down list, if the IT guys are worth their salt) to save on time. The geologist must also 'log' what samples are to be sent to the lab for assay analysis. This is what the geologist thinks either looks like ore or contains features that might contain significant grades (Quartz veins or sulphides for a gold deposit for example). Many management crews will insist on the rock codes aligning with assay sample intervals however this is not always the case.

I am trying to assign an ore value (Au) to a rock code. In two separate files, first one is a list of ore assay samples from drillholes at approximate 1 or 2 m intervals. The intervals are not always 1 or 2m but in general it is attempted to keep the samples at similar intervals for compositing purposes. The second file is the lithology rockcodes from the same drillholes. Depending on the width of lithology determines the interval width.

Im trying to get the average assay grade for each lithology (rockcode) interval. This would mean the sum of all the assay values over the width of the lithology interval. The complicated part is when the assay samples overlap across a lithology sample boundary, meaning I have to incorporate into the code the ratio of the first assay sample that is divided by the lithology boundary. I guess i can see how to do it, but cannot get the coding to do what I want.

So for example, you may have an interval of granite that is 8m wide, and say it contained 4 samples of 2m width and they conveniently aligned with the beginning and end of the granite interval, I could simply sum these 4 Au values together and say this is the grade of the granite over 8m thickness. However, If the first assay sample started at .5m before the initial granite interval, I would have to weight the assay sample to 1.5m (multiply the assay value by .75) and this would be the first assay of the granite. Followed by 2 whole 2m samples with the final assay being weighted by .25 (because only .5m of this sample sits within the granite)

The fraction of the interval above the overlap must be weighted against the ratio of the thickness in order to figure out what Au value it is and then added to the Au value that is either enveloped within the width of the rockcode interval or added to its lower fraction if there is no assay interval that is completely within the rockcode interval.

Attached are shortened versions of my files with just three holes included to make it a little easier. The drillholes are CSH-11-001, CSH-11-002, CSH-11-003.

You can see from the section titled 'do the smart thing' its a coding problem.

Here are both files:

Here is a screenshot of both files with the first 19.5m done manually in the centre..

[url src="" [/u][/url]

Here is my code so far....

Here is a screenshot of the first 20m completed manually. I hope its clear what im trying to do, if not i can put together a graphic showing how I calculated the rockcode grade intervals manually...


'Declare variables that represent each column in the .csv file

nHoleID = 1
nFrom = 2
nTo = 3
nLength = 4
nCode = 5
nGrade = 6

nHoleID = 1
nFrom = 2
nTo = 3
nLength = 4
nAu = 5

'count the number of records in the file
LithCount = 0
OPEN "export_lithos.csv" for INPUT as #count
   WHILE EOF(#count) = 0
        LithCount = LithCount + 1
        line input #count, n$
CLOSE #count
PRINT "There are ";LithCount;" items of data in this file."
print "-------------------------------------------"
LithCount = LithCount - 1

'count the number of records in the file #2
AssCount = 0
OPEN "export_assay.csv" for INPUT as #count1
   WHILE EOF(#count1) = 0
        AssCount = AssCount + 1
        line input #count1, n$
CLOSE #count1
PRINT "There are ";AssCount;" items of data in this file."
print "-------------------------------------------"
AssCount = AssCount - 1

'Dimension the arrays to hold the data

dim LithHeader$(6)
dim LithHoleID$(LithCount)
dim LithData(LithCount, 6)

dim AssHeader$(5)
dim AssHoleID$(AssCount)
dim AssData(AssCount, 5)

'Open file and input data header and data from columns into array(s)

OPEN "export_lithos.csv" for INPUT as #readme

line input #readme, Header$

for i = 1 to LithCount
    line input #readme, data$
    LithHoleID$(i) = word$(data$,1, ",")
    LithData(i,nFrom) = val(word$(data$,2, ","))
    LithData(i,nTo) = val(word$(data$,3, ","))
    LithData(i,nLength) = val(word$(data$,4, ","))
    LithData(i,nCode) = val(word$(data$,5, ","))
    LithData(i,nGrade) = 0


close #readme

'Open file and input data header and data from columns into array(s)

OPEN "export_assay.csv" for INPUT as #readthis

line input #readthis, HeaderAss$

for i = 1 to AssCount
    line input #readthis, data$
    AssHoleID$(i) = word$(data$,1, ",")
    AssData(i,nFrom) = val(word$(data$,2, ","))
    AssData(i,nTo) = val(word$(data$,3, ","))
    AssData(i,nLength) = val(word$(data$,4, ","))
    AssData(i,nAu) = val(word$(data$,5, ","))


close #readthis

'Create a new output file

open "Lithos_grade_result.csv" for output as #outfile

'Print the header with new columns name(s) to the output file

print #outfile, Header$;","; "nGrade"

'do the smart thing

while i < LithCount
    while LithHoleID$(i) = AssHoleID$(j)
        if LithData(i,nFrom) = AssData(j,nFrom) and LithData(i,nTo) = AssData(j,nTo) then
           LithData(i,nGrade) = AssData(j,nAu)
        end if
        j = j + 1
        i = i + 1
            if LithData(i,nFrom) = AssData(j,nFrom) and LithData(i,nTo) < AssData(j,nTo) then
                LithData(i,nGrade) = (LithData(i,nLength) / AssData(j,nLength)) * AssData(j,nAu)
            end if

            j = j + 1

                if AssData(j,nFrom) > LithData(i,nTo) then
                    fractionAbv = (AssData(j,nFrom) - LithData(i,nTo)) * AssData(j-1,nAu)
                        while AssData(j,nTo) < LithData(i,nTo)
                            Avelength = Avelength + AssData(j,nLength)
                            totalore = totalore + AssData(j,nAu)
                            j = j + 1
                        fractionAbv = (AssData(j,nTo) - LithData(i,nTo)) * AssData(j,nAu)
                        LithData(i,nGrade) = (totalore/Avelength) + fractionAbv + fractionBlw
                end if

i = i + 1

fractionAbv = 0
fractionBlw = 0
totalore = 0
Avelength = 0


'Print all data from arrays to outfile appending calculations to latter columns

for i = 1 to LithCount

    print #outfile, LithHoleID$(i);",";
    print #outfile, LithData(i,nFrom);",";
    print #outfile, LithData(i,nTo);",";
    print #outfile, LithData(i,nLength);",";
    print #outfile, LithData(i,nCode);",";
    print #outfile, LithData(i,nGrade)


close #outfile



Re: simple math but confusing coding
Post by blessed2b3rd on Oct 20th, 2015, 9:38pm

I realize that using percentages may complicate things. If one was to use the weighted grades by multiplying by the respective lengths and then divide by the sum of the length to get the average grade it might make the coding a little easier
Re: simple math but confusing coding
Post by tenochtitlanuk on Oct 21st, 2015, 05:18am

Struggling a bit to understand the physics of this. Wrote a quick graphic prog. to clarify for me and others. Shows the rock layers from the litho file and the assay values. Use 'Show Image' if your browser shows it scaled down.

User Image

I'd have expected to want the layer's vertical extent to be multiplied by the assay value within it. But if the assay value is on a boundary, is it supposed to represent the rock above or below? If it is to represent the rock BELOW then in the case of 'start and finish on a boundary' the last value should be ignored as belonging in fact to the next stratum down???

PS it also briefly held me up that the rock type classifications are all saved with a leading space for some reason!

PPS It is indeed interesting to see real-world applications of this kind of data-mining.

PPPS. LB 4.5 makes reading csv files even easier-
-Added an INPUTCSV statement for reading comma delimited string values
 from a file with sensitivity to quoted values.

  open "someFile.txt" for input as #1
  inputcsv #1, a$, b$, c$, d$
  print a$
  print b$
  print c$
  print d$
  close #1 

Code used to produce strata graphic..

WindowWidth  = 520
WindowHeight = 700

open "Miner49er" for graphics_nsb as #wg

global assayOpen, lithosOpen

#wg "goto 0 20 ; down ;fill 200 180 170 ; goto 700 20 ; trapclose quit ; flush"

oldLevel = 0
vScale   = 9

open "export_lithos.csv" for input as #lithos
   lithosOpen =1
   while eof( #lithos) =0 and Level <700 /vScale
        line input #lithos, n$
        if word$( n$, 1, ",") ="CSH-11-001" then
            Level =int( 0.5 +val( word$( n$, 2, ",")))
            rockCode$ =word$( n$, 5, ",")
            col$ =lookUp$( rockCode$)
            print Level, rockCode$, col$
            #wg "backcolor "; col$
            #wg "place 20 "; 20 +vScale *oldLevel; " ; boxfilled 50 "; 20 +vScale *Level
            'print "#wg 'place 20 "; 20 +vScale *oldLevel; " ; boxfilled 50 "; 20 +vScale *Level; "'"
            oldLevel =Level
        end if

close #lithos
lithosOpen =0

print: print

Level =0

#wg "color red ; size 4 ; flush ; discard ; down"

open "export_assay.csv" for input as #assay
   assayOpen =1
   while eof( #assay) =0 and Level <700 /vScale
        line input #assay, n$
        if word$( n$, 1, ",") ="CSH-11-001" then
            Level =int( 0.5 +val( word$( n$, 2, ",")))
            Au =val( word$( n$, 5, ","))
            'print "#wg 'set 40 "; 20 +vScale *Level; "'"
            #wg "up ; goto 50 "; 20 +vScale *Level; " ; down ; goto "; 50 +Au *100; "  "; 20 +vScale *Level
        end if
close #assay
assayOpen =0

#wg "flush ; getbmp scr 0 0 519 719"
bmpsave "scr", "Hole CSH-11-001.bmp"

wait    '   _______________________________________________________

function lookUp$( r$)   '   so can colour-code graphic output....
    select case r$
        case " S"
            lookUp$ ="200 200 200"
        case " mC2"
            lookUp$ ="200 100 0"
        case " mC1"
            lookUp$ ="100 200 50"
        case " mAC"
            lookUp$ ="70 120 30"
        case " mC3"
            lookUp$ ="100 0 50"
        'case " mC1"
            'lookUp$ ="50 0 250"
        case else
            lookUp$ ="120 120 190"
    end select
end function

sub quit h$
    close #wg
    if lithosOpen then close #lithos
    if assayOpen  then close #assay
end sub

Re: simple math but confusing coding
Post by blessed2b3rd on Oct 21st, 2015, 10:37am

Thanks Teno,

Good eye with the leading spaces, not sure how that error slipped through.

Ive adapted your graphics code to 'un-round' the interval boundaries to make it more accurate. This way one isn't fooled into thinking that the boundaries are nicely overlapping one another. See my picture to show how to calculate one of the intervals.

Ignore the screenshot in my previous post as I calculated incorrectly.

User Image
Re: simple math but confusing coding
Post by tenochtitlanuk on Oct 22nd, 2015, 08:08am

I'm still mixed up. I'd expect the result to be calculated with the assay values of 0.058 and 0.845 rather than 0.086 and 0.058 huh
User Image
PS away from base for a few days...
Re: simple math but confusing coding
Post by blessed2b3rd on Oct 22nd, 2015, 09:28am

OMG, yea you are right, I pushed the first assay value (0-0.4m) onto the second interval